AFIT/GNE/ENP/93M-6 


•S  i 


CONSTRUCTION  AND  TESTING  OF  A 
NEUTRON  AND  GAMMA  SPECTROMETRY 
SYSTEM  USING  PULSE  SHAPE  DISCRIMINATION 
WITH  AN  ORGANIC  SCINTILLATOR 

THESIS 

Robert  S.  Pope,  Second  Lieutenant,  USAF 


OTIC 

ELECTE 
APR  05  1993 

E 


AFIT/GNE/ENP/93M-6  - 


Reproduced  From 
Best  Available  Copy 


4  02  057 


Approved  for  public  release,  distribution  unlimited 

^000 W  0% 


AFIT/GNE/ENP/93M-6 


CONSTRUCTION  AND  TESTING  OF  A  NEUTRON  AND  GAMMA 
SPECTROMETRY  SYSTEM  USING  PULSE  SHAPE  DISCRIMINATION 
WITH  AN  ORGANIC  SCINTILLATOR 


THESIS 


Presented  to  the  Faculty  of  the  School  of  Engineering 


of  the  Air  Force  Institute  of  Technology 


Air  University 
In  Partial  Fulfillment  of  the 
Requirements  for  the  Degree  of 
Master  of  Science 


Robert  S.  Pope,  B.S. 
Second  Lieutenant,  USAF 


Accesion  For 

NTIS  CRA&I 
OTIC  TAB 
Unannounced 
Justification 

□ 

By . 

Distribution  / 

. . . — 

Availability  Codes 


Dist 


Avail  and/or 
Special 


/ 


March  1993 


Approved  for  public  release;  distribution  unlimited 


Acknowledgments 


The  contributions  of  many  people  have  made  this  thesis  possible.  I  would  tirst  like 
to  thank  my  advisor,  LtCol  (Dr.)  Richard  Hartley,  and  committee  members.  Dr.  Nolan 
Hertel  and  Dr.  George  John.  LtCoI  Hartley  provided  a  great  deal  of  time  and  guidance 
throughout  the  course  of  this  project,  which  began  with  my  Research  Apprenticeship  work 
more  than  eight  months  ago.  Dr.  Hertel  offered  guidance  with  the  PSD  setup  of  the 
spectrometry  system  and  provided  many  helpful  suggestions  on  the  modifications  of 
FORIST  and  ZSHIFT.  Dr.  John,  who  taught  the  Nuclear  Physics  and  Nucleonics 
Instrumentation  courses  at  the  Air  Force  Institute  of  Technology,  provided  me  with  the 
knowledge  base  necessary  to  begin  this  research  project. 

I  would  also  like  to  thank  Doug  Burkholtz  and  Jack  Phillips  of  the  AFIT  small 
computers  section  for  all  of  their  help  with  the  VAX  VMS  mainframe.  Additionally,  I 
would  like  to  thank  my  classmate,  Capt  Sean  Miller,  who  worked  with  me  during  the 
Research  Apprenticeship  project  and  provided  his  encouragement  and  a  listening  ear 
throughout  the  thesis  work.  I  would  like  to  especially  thank  Bob  Hendricks,  the 
laboratory  supervisor  at  the  Air  Force  Nuclear  Engineering  Center,  who  provided 
invaluable  support  and  assistance  in  the  laboratory. 

Finally,  there  are  two  literature  sources  which  deserve  special  mention.  I  highly 
recommend  the  book  Radiation  Detection  and  Measurement  by  Glenn  F.  Knoll  to  anyone 
who  plans  to  do  an  experimental  thesis  in  the  field  of  Nuclear  Physics.  While  the  literature 
search  performed  for  this  thesis  was  extensive,  nearly  ail  inquires  began  with  Knoll. 
Additionally,  the  best  periodical  source  of  experimentally-oriented  literature  that  was 
found  during  the  course  of  this  project  was  the  journal  Nuclear  Instruments  and  Methods. 


Page 


Table  of  Contents 


Acknowledgments  .  ii 

List  of  Figures  .  v 

List  of  Tables  . viii 

Abstract  . . .  ix 

I.  Introduction  .  1 

Background  .  1 

Scope  .  3 

Sequence  of  Presentation  .  4 

II.  Theory  . 5 

Pulse  Shape  Discrimination  (PSD)  .  5 

Spectrum  Unfolding  .  10 

Energy  Calibration  .  22 

III.  Description  of  Apparatus  .  26 

Detector  . 26 

Electronics  . 30 

Shadow  Shield  .  34 

Radiation  Sources  .  35 

Computer  Software  .  37 

IV.  Experimental  Procedure  .  38 

Laboratory  Arrangement  .  38 

Electronics:  Set-up  and  Optimization  .  40 

Data  Collection  Procedure  .  52 

Data  Processing  and  Unfolding  Procedure  .  66 


V.  Experimental  Results  . 

Test  of  Modified  FORIST  Code 

PuBe  Source  #T022  . 

PuBe  Source  #M- 1170  . 

Na-22  Source  #T028  . 

VI.  Summary  and  Conclusions  . 

Appendix  A:  Equipment  List  . 

Appendix  B;  The  Plutonium-Beryllium  Source 
Appendix  C:  CALIBER,  TK  Computer  Model 
Appendix  D:  ZSHIFT.FOR  Computer  Code  ... 
Appendix  E;  FORIST.FOR  Computer  Code  .  . 
Appendix  F:  GRAPHIC S.TK  Computer  Model 

Appendix  G:  ZSHIFT.OUT  Data  File  . 

Appendix  H;  SPECTRA.  IN  Data  File  . 

Appendix  I:  FORIST.OUT  Data  File  . 

Appendix  J:  TOGRAPH.OUT  Data  File  . 

Bibliography  . 

Vita  . 


Page 

72 

72 

77 

80 

82 

85 

86 
92 
94 
98 

122 

123 

133 

137 

179 

193 

196 


List  of  Figures 


Figure  Page 

1 .  Appro:iimate  energy  ranges  of  several  neutron  spectrometers .  2 

2.  Energy  levels  of  an  organic  molecule  with  .t-electron  structure .  7 

3.  Time  dependence  of  scintillation  pulses  to  radiation  of  different  types .  9 

4.  Ideal  energy  distribution  of  recoil  protons  .  11 

5.  Distortion  to  the  recoil  proton  response  function  due  to 

scattering  with  carbon .  13 

6.  Distortion  to  the  recoil  proton  response  function  due  to 

finite  detector  resolution  .  13 

7.  (a)  Pulse  height  vs  neutron  energy;  (b)  Derivative  of  (a); 

(c)  Distortion  of  rectangular  proton  recoil  spectrum  due  to 

nonlinear  response  of  scintillator .  14 

8.  Distortion  to  the  recoil  proton  response  function  due  to 

multiple  scatters  with  hydrogen .  16 

9.  Ideal  response  curves  of  neutrons  of  several  discrete  energies .  17 

10.  Total  ideal  response  function  of  neutrons  of  several  discrete  energies .  18 

1 1 .  Energy  distribution  of  Compton  recoil  electrons .  22 

12.  22Na  Compton  recoil  electron  spectrum  from  literature .  23 

13.  Diagram  of  detector  used  in  this  research  project  .  27 

14.  Block  diagram  of  electronic  equipment  for  spectrometry  system  .  30 

15.  Laboratory  arrangement  of  source,  detector,  and  shadow  shield .  39 

16.  Electronics  configuration  for  PSD  work .  40 


1 7.  Experimentally-determined  optimum  PSD  for  this  spectrometry  system 


41 


Page 


Figure 

18.  Oscilloscope  trace  of  output  from  the  preamp  showing 

several  signals  from  the  PuBe  #T022  y/n  source .  43 

19.  Oscilloscope  trace  of  output  from  preamp  showing  single  signal  from 

PuBe  #T022  with  magnified  horizontal  scale  to  show  detail  of  pulse  .  44 

20.  Oscilloscope  trace  of  bipolar  output  from  DLA  showing 

single  signal  from  PuBe  #T022  . .  45 

2 1 .  Oscilloscope  trace  of  linear  output  from  PSA  showing 

a  single  signal  from  PuBe  #T022  .  46 

22.  Oscilloscope  trace  of  linear  output  from  PSA  using  5x  compression 
feature  on  scope  to  show  4  gamma  pulses  (-4.3  V) 

and  1  neutron  pulse  (-7.6  V) .  47 

23.  Oscilloscope  trace  showing  the  inputs  to  the  ADC  used  for  setting  PSD .  49 

24.  Oscilloscope  trace  of  inputs  to  the  ADC  for  recording  energy  spectra .  51 

25.  Laboratory  measurement  of  #T028  ^2^3  Compton  electron  spectrum 

at  high  gain  setting .  55 

26.  Laboratory  measurement  of  #T028  22Na  Compton  electron  spectrum 

at  low  gain  setting .  56 

27.  Three  tables  created  by  CALIBER. TK .  58 

28.  Plot  of  peak  positions  vs.  amplifier  gain  to  verify  linearity 

of  the  data  in  Figure  27 .  59 

29.  Plot  of  gamma  1  peak  and  half-height  values  vs.  amplifier  gain 

to  verify  linearity .  60 

30.  Part  of  CALIBER.TK  variable  sheet  showing  gain  and 

zero-intercept  results .  61 

3 1 .  PuBe  #M- 1170  Compton  recoil  electron  high  gain  foreground  spectrum .  63 

32.  PuBe  #M- 1 1 70  Compton  recoil  electron  high  gain  background  spectrum  ....  64 


VI 


Figure  Page 

33.  PuBe  #M- 1 1 70  Compton  recoil  electron  low  gain  foreground  spectrum .  65 

34.  PuBe  #M-1 170  Compton  recoil  electron  low  gain  background  spectrum .  66 

35.  DT  neutron  spectn.'m  used  to  verify  correct  operation 

of  modified  FORIST  code . .  73 

36.  Experimentally-determined  PuBe  #T022  neutron  spectrum .  75 

37.  Experimentally-determined  PuBe  #T022  gamma  spectrum  .  76 

38.  Experimentally-determined  PuBe  #M-1 170  neutron  spectrum  .  78 

39.  Experimentally-determined  PuBe  #M-1 170  gamma  spectrum .  79 

I 

40.  Experimentally-determined  22Na  #T028  gamma  spectrum .  81 

41.  Plutonium  decay  series .  87 

42.  Theoretical  spectrum  of  a  239puBe  source  showing  total  and  partial 

neutron  spectra  corresponding  to  continuum,  ground,  first,  and  second 
excited  states  of  *2c  .  91 


vii 


List  of  Tables 


Table  Page 

1.  Comparison  of  NE-2 13  and  BC-501 A  scintillators  .  28 

2.  Comparison  of  RCA  8610-A  and  Hamamatsu  R329-02 

photomultiplier  tubes  .  29 

3.  Summary  of  electronics  settings  .  53 

4.  FORISl  parameters  for  neutron  and  gamma  unfolding .  70 

5.  Inputs  to  ZSHIFT  for  FORIST  verification  run  . 73 

6.  Inputs  to  ZSHIFT  for  PuBe  Ul  022  neutron  run  .  74 

7.  Comparison  of  true  PuBe  neutron  energies  to  measured 

PuBe  #T022  spectrum . 75 

8.  Inputs  to  ZSHIFT  for  PuBe  #T022  gamma  run  .  76 

9.  Inputs  to  ZSHIFT  for  PuBe  #M-1 170  neutron  run  .  77 

10.  Comparison  of  true  PuBe  neutron  energies  to  measured 

PuBe  #M-1 170  spectrum  . 78 

1 1 .  Inputs  to  ZSHIFT  for  PuBe  #M- 1 1 70  gamma  run .  80 

12.  Inputs  to  ZSHIFT  for  ^^Na  #T028  gamma  run . .  80 

13.  Equipment  list .  85 

14.  Description  of  contents  of  FORIST.OUT  file .  137 


AFTT/GNE/ENP/93M-6 


Abstract 


The  goal  of  this  thesis  was  to  construct  and  test  a  neutron  detection  system  to 
measure  the  energy  spectrum  of  neutrons  in  the  presence  of  gammas  for  neutrons  in  an 
energy  range  from  1  to  14  MeV.  To  accomplish  this,  a  neutron  spectiometer  based  on  the 
process  of  pulse  shape  discrimination  (PSD)  was  constructed,  in  which  a  xylene-based 
scintillator,  NE-213,  was  used  as  the  detector  material.  The  primary  neutron/gamma 
sources  used  were  78-mCi  and  4.7-Ci  ^^^PuBe  sources,  while  4.7-pCi  and  97.6-pCi  22Na 
gamma  sources  were  used  for  energy  calibration  and  additional  testing  of  the  detector. 
Proton  recoil  spectra  and  Compton  electron  recoil  spectra  were  unfolded  with  the  neutron 
and  gamma  unfolding  code  FORIST  to  generate  the  incident  neutron  and  gamma  spectra, 
respectively.  FORIST,  which  was  written  for  a  CDC  computer,  was  modified  to  run  on  a 
VAX  6420. 

The  experimental  spectra  were  compared  to  those  in  the  scientific  literature.  It 
was  found  that  the  locations  of  the  peaks  in  the  239puBe  spectrum  agreed  with  the 
literature  to  within  8.3%,  the  239puBe  gamma  spectrum  agreed  to  within  0.7%,  while  the 
22Na  gamma  spectrum  agreed  exactly.  Uncertainties  in  the  detection  system  and 
unfolding  procedure  are  on  the  order  of  5  to  10%,  so  these  results  are  reasonable. 

This  thesis  was  written  as  a  summary  of  the  relevant  literature  and  as  a  user's  guide 
to  the  neutron  spectrometer  developed  in  this  project.  To  this  end,  this  paper  first 
develops  the  theory  of  PSD  techniques,  computer  unfolding,  and  energy  calibration;  then 
describes  the  hardware,  software,  and  experimental  procedures;  and  in  the  end  presents 
the  experimental  results  and  compares  them  to  the  literature  results  as  validation  of  the 
system. 


CONSTRUCTION  AND  TESTING  OF  A  NTUTRON  GAMMA 
SPECTROMETRY  SYSTEM  USING  PULSE  SHAPE  DISCRIMINATION 
WITH  AN  ORGANIC  SCINTILLATOR 


I.  Introduction 

Background 

Neutron  detection  and  spectrometry  are  very  important  to  nuclear  research,  but 
because  the  neutron  is  an  uncharged  particle,  detecting  neutrons  and  measuring  their 
energy  is  not  a  trivial  task.  Because  the  energy  range  of  neutrons  can  cover  !4  decades 
(lO"^  to  10*  eV),  many  diflferent  neutron  detection  and  spectrometry  methods  have  been 
developed,  most  of  which  only  cover  a  portion  of  ihe  full  neutron  energy  range. 

Nearly  all  present  neutron  spectrometry  methods  are  based  on  one  of  four  general 
principles; 

1 .  measurement  of  the  energies  of  charged  particles  produced  by  neutron 
interactions,  e.g.,  nuclear  emulsions,  damage  track  detectors,  recoil  proton 
telescopes,  proportional  counters,  sandwich  spectrometers,  lithium  glass 
scintillators,  and  organic  scintillators; 

2.  measurement  of  the  flight  time  of  neutrons  over  a  given  distance,  e.g.,  TOF 
spectrometers; 

3.  measurement  of  attenuation  and  moderation  of  neutrons  in  various 


thicknesses  of  hydrogenous  materia!,  e.g.,  Bonner  spheres;  and 
4.  measurement  of  activation  or  fission  reactions  in  a  set  of  detectors  having 
different  energy  responses,  e.g..,  threshold  detectors. 

The  approximate  energy  ranges  covered  by  these  diffeient  spectrometer  types  a^-e  shown 
in  Figure  1  [8:95-96], 
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Figure  1.  Approximate  energy  ranges  of  several  neutron  spectrometers  [8:96] 


The  detector  constructed  for  this  research  project  operates  on  the  first  of  the  four 
principles:  the  measurement  of  the  energies  of  charged  particle^  produced  by  neutron  and 

j 

gamma  interactions.  The  detection  material  is  a  xylene-based  liquid  organic  scintillator, 
often  called  NE-213  in  the  literature.  This  organic  scintillator  dejtects  both  incoming 
neutrons  and  gammas.  Neutron-gamma  discrimination  is  accomplished  by  the  method  of 
piilse  shape  discrimination  (PSD),  which  uses  a  timing  circuit  that  operates  on  the 
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difference  in  the  decay  component  of  the  output  signal  (the  shape  of  the  pulse)  from 
neutrons  and  gammas.  The  timing  information  is  then  used  to  gate  the  energy 
information,  so  that  the  user  can  select  whether  the  neutron  or  gamma  energy  is  recorded, 
PSD  is  useful  over  the  entire  l-to-H-ivIeV  range  of  interest.  Pulse  shape  discrimination  is 
explained  more  fully  in  the  "Theory"  section. 

The  organic  scintillator  does  not  record  the  neutron  and  gamma  energy  spectra 
directly.  Rather,  incoming  neutrons  interact  with  the  hydrogen  nuclei  in  the  scintillator  to 
produce  recoil  protons,  v/hile  the  incoming  gammas  interact  with  atomic  electrons  to 
produce  Compton  recoil  electrons.  Re-creating  the  incoming  neutron  and  gamma  spectra 
requires  unfolding  the  recoil  spectra  recorded  by  the  detector.  Unfolding  is  a 
mathematical  operation  in  which  the  recoil  information  is  placed  into  a  vector  and 
multiplied  with  an  appropriate  response  matrix  to  yield  a  vector  which  contains  the 
calculated  incident  n^’.utron  or  gamma  spectrum.  Spectrum  unfolding  is  discussed  in  detail 
in  the  "Theory"  section.  The  FORIST  (FERDOR  with  Optimized  Resolution  using  an 
Iterative  Smoothing  Technique)  [26;-]  unfolding  code,  modified  to  run  on  the  VAX  6420, 
is  used  to  unfold  the  experimental  recoil  spectra  recorded  with  the  detector. 


Scope 


The  type  of  detection  system  used  in  this  project  is  not  new.  Verbinski  et  al. 
[29:-;6:-]  did  the  definitive  work  with  an  NE-213  neutron  spectrometer  in  1968.  In  that 
work,  Verbinski's  team  used  monoenergetic  neutron  sources  to  construct  the  neutron 
response  matrix  for  their  detector.  That  response  matrix  forms  the  basis  for  the  neutron 
response  matrix  used  with  FORIST.  They  also  established  protocols  in  equipment  and 
procedure  that  many  others  subsequently  followed  so  that  Verbinski's  response  matrix 
would  apply  to  their  work.  To  be  able  to  use  the  FORIST  unfolding  code,  the  detector 
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built  for  the  Air  Force  Nuclear  Engineering  Center  during  this  thesis  project  also  follows 
Verbinski's  protocols.  A  comparison  of  the  two  detectors  is  made  in  the  "Description  of 
Apparatus"  section. 

The  focus  of  this  thesis  research  is,  then,  not  to  develop  a  new  kind  of  neutron 
spectrometer,  but  to  construct  and  test  an  accepted  spectrometer  design  and  then  to 
provide  a  complete  user’s  manual  for  the  system.  To  that  end,  this  thesis  completely 
describes  the  apparatus  and  computer  codes  used  and  then  gives  a  complete  procedure  to 
follow  in  order  to  construct  and  use  this  system.  A  search  of  the  literature  indicated  that, 
while  this  kind  of  spectrometer  has  been  in  use  for  more  than  twenty  years,  there  was  no 
complete  guide  for  the  system  geared  to  the  new  user.  The  intent  of  this  thesis  is  to 
provide  such  a  manual. 

Sequence  of  Presentation 

The  rest  of  this  thesis  is  divided  into  five  major  sections.  The  first  section 
discusses  the  general  theovy  behind  pulse  shape  discrimination  and  spectrum  unfolding. 
The  second  section  describes  the  apparatus  and  computer  software  used  during  this 
project.  The  third  section  presents  a  step-by-step  guide  to  assembling  and  using  the 
spectrometer  constructed  for  this  research,  from  initial  setup  through  data  collection  and 
processing.  The  fourth  section  presents  the  experimental  results  from  the  project.  These 
results  include  verification  of  the  correct  operation  of  the  modified  FORIST  code  and 
validation  of  the  spectrometry  sytstem.  This  validation  was  performed  by  comparing 
exprimental  spectra  from  239puBe  and  22Na  to  spectra  from  the  literature.  The  fifth  and 
final  section  summarizes  the  project  and  presents  recommendations  for  future  work  with 
this  detection  system. 
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11.  Theory 


Pulse  Shape  Discrimination  fPSDl 

Organic  scintillation  detectors  have  been  in  use  for  more  than  forty  years  [4:477], 
These  detectors  respond  both  to  incident  gammas  and  neutrons.  Distinguishing  between 
neutrons  and  gammas  is  the  function  of  pulse  shape  discrimination  (PSD).  In  order  to 
explain  PSD,  one  must  first  understand  the  response  of  organic  scintillators,  such  as 
NE-213,  to  incident  radiation.  To  that  end,  this  section  first  describes  how  the  uncharged 
neutrons  and  gammas  are  detected  within  the  scintillator.  It  then  presents  a  discussion  of 
the  scintillation  mechanism  in  organic  scintillation  detectors,  which  is  used  to  explain  the 
theory  of  PSD. 

Scintillation  is  the  process  in  which  incident  ionizing  radiation  upon  a  material 
produces  the  emission  of  light.  There  are  two  types  of  scintillators:  organic  and  inorganic. 
In  general,  inorganic  scintillators  are  favored  for  gamma-ray  spectroscopy  because  of  their 
high-Z  constituents.  Though  comparatively  less  efficient,  organic  scintillators  are  selected 
for  fast  neutron  detection  because  of  advantages  in  energy  resolution  and  pulse  shape 
discrimination.  Organic  scintillators  are  composed  of  compounds  of  hydrogen  and  carbon 
and  are  available  in  crystalline,  plastic,  or  liquid  form.  The  scintillator  used  in  this  project 
is  a  binary  liquid  scintillator  in  which  a  small  quantity  of  fluor  is  dissolved  in  a  large 
quantity  of  organic  solvent  (xylene).  The  organic  solvent  participates  in  the  reactions  with 
incident  radiation  and  then  transfers  the  resulting  energy  to  the  fluor,  which  then  de- 
excites  by  emitting  light  [4:479]. 

Neutrons  entering  the  scintillator  interact  by  elastically  scattering  with  hydrogen 
and  carbon  nuclei.  This  elastic  scattering  reaction  transfers  some  or  all  of  the  incident 


5 


neutron's  energy  to  the  target  nucleus,  which  then  recoils.  At  neutron  energies  above 
about  8  MeV,  -2C((X,n)  reactions  are  also  possible,  which  create  multiple  recoiling 
particles.  The  recoil  nucleus  or  I2c(a,n)  reaction  products  then  collide  with  the  fluor, 
transferring  their  excitation  energy.  The  fluor  subsequently  de-excites  by  emitting  light 
[4:485-486,17:530-540]. 

Gamma  photons  entering  the  scintillator  interact,  not  with  the  scintillator  nuclei, 
but  with  the  atomic  electrons.  Due  to  the  low-Z  scintillator  constituents,  they  interact 
almost  entirely  by  Compton  scatter.  The  Compton  interaction  creates  a  Compton  recoil 
electron  [17:54]  that  then  transfers  its  energy  to  the  fluor,  which  subsequently  de-excites 
through  fluorescence. 

Organic  scintillators  have  a  molecular  structure  known  as  n-electron  structure. 
This  structure  consists  of  planar  molecules  built  up  mainly  from  benzene  rings.  A  key 
feature  of  the  molecular  structure  is  the  presence  of  extended  groupings  of  conjugate"^ 
double  bonds,  based  largely  on  linkages  between  unsaturated  carbon  atoms.  This  nor 
saturation  implies  that  only  two  or  three  of  the  four  valence  electrons  of  carbon  are 
strongly  localized  within  the  molecular  structure.  The  remaining  valance  electrons  art 
delocalized  within  the  molecule  and  are  not  associated  with  any  one  particular  atom, 
occupy  7t-molecular  orbitals  which  extend  over  the  conjugated  region.  The  7t-electro 
states  are  of  particular  interest  because  transitions  between  these  states  lead  to  the 
luminescence  observed  in  the  scintillation  process  [4:479],  The  7i-electronic  energy  It  . 
of  such  a  molecule  are  illustrated  in  Figure  2. 

Scintillation  occurs  in  a  molecule  with  ;r-electron  structure  in  the  following 
manner.  In  a  molecule  at  thermal  equilibrium,  the  electrons  are  in  the  lowest  vibratio 
state  of  the  ground  electronic  state,  Sqo.  Electrons  are  excited  into  a  higher  energy  Ir 
when  energy  is  absorbed  by  the  molecule.  They  can  be  excited  into  one  of  the  single, 
states  (spin  0)  or  into  one  of  the  triplet  states  (spin  1).  The  first  subscript  denotes  as 
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Singlet 


Triplet 


Figure  2.  Energy  levels  of  an  organic  molecule  with  n-electron  structure  [17:217] 
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quantum  number  and  the  second  a  vibrational  quantum  number.  The  energy  of  the  Sj 
state  is  approximately  3-4  eV  above  Sq  and  the  energy  difference  between  higher  states  is 
generally  smaller.  The  energy  differences  between  vibrational  states,  Sj^,  are  on  the  order 
of  0. 1 5  eV  apart.  Because  this  is  relatively  large  compared  to  thermal  temperatures 
(0.025  eV),  nearly  all  molecules  are  in  the  Sqq  state  at  room  temperature.  [17:216-217] 
When  an  energetic  particle  such  as  a  recoil  proton  or  Compton  recoil  electron 
collides  with  the  fluor  molecule,  some  of  the  kinetic  energy  is  absorbed  by  the  molecule. 
This  is  represented  by  the  up  arrows.  Electrons  which  are  excited  into  the  higher  singlet 
states  quickly  de-excite  to  the  Sjq  level  through  radiationless  internal  conversions.  Any 
electrons  in  Sj  above  Sjq  are  not  in  thermal  equilibrium  with  their  neighbors  and  also 
quickly  lose  the  excess  vibrational  energy.  The  net  effect  of  a  single  energetic  particle 
passing  through  the  material  is  to  produce  a  large  number  of  molecules  with  excited 
electrons  all  in  the  same  energy'  state.  The  number  of  electrons  in  the  higher  energy  state 
and  the  consequent  fluorescence  are  proportional  to  the  kinetic  energy  of  the  recoiling 
particle  that  passed  through  scintillating  material  [17:217].  The  scintillating  material  is 
said  to  fluoresce  when  the  electrons  in  the  Sjq  state  de-excite  to  the  ground  state  or  one  of 
the  ground  vibrational  states.  These  prompt  decays  are  indicated  by  the  down  arrows  in 
Figure  2.  The  time  it  takes  for  a  material  to  fluoresce  once  it  has  been  excited  is  on  the 
order  of  a  few  nanoseconds  [17:226], 

Electrons  which  are  excited  through  inter-system  crossings  to  a  triplet  state  must 
de-excite  in  a  different  manner.  When  electrons  in  the  Tj  state  de-excite  to  the  ground 
state  or  one  of  the  ground  vibrational  states,  the  light  emitted  in  this  delayed  decay  is 
called  phosphorescence.  However,  a  decay  from  the  T  j  triplet  state  directly  to  the  Sqo 
ground  state  is  not  a  quantum-mechanically-preferred  mode  of  decay  and  hence 
contributes  little  to  the  light  output  [4:480],  The  preferred  mode  of  de-excitation  for  a 
molecule  in  a  triplet  state  is  to  participate  in  a  bimolecular  interaction  with  another 


molecule  in  the  triplet  state.  Such  an  interaction  places  one  molecule  in  the  first  singlet 
state  (Sj)  and  the  other  in  the  ground  state  (Sq)  The  molecule  in  the  Sj  state  can  then  de- 
excite  in  the  normal  way,  leading  to  delayed  fluorescence,  which  takes  on  the  order  of  a 
few  hundred  nanoseconds  [17:227], 

The  majority  of  the  light  yield  of  the  scintillator  occurs  in  the  prompt  fluorescence 
component,  but  the  delayed  fluorescence  component  is  significant  because  the  fraction  of 
light  that  appears  in  the  delayed  component  depends  on  the  nature  of  the  exciting  particle. 
The  delayed  fluorescence,  then,  provides  the  basis  for  pulse  shape  discrimination  [17:226], 
In  general,  as  the  mass  of  the  incident  particle  increases,  the  number  of  fluor  molecules  de¬ 
exciting  through  delayed  fluorescence  also  increases  and  hence  increases  the  delayed 
component  of  the  tesulting  pulse.  This  is  shown  in  Figure  3.  Pulse  shape  discrimination, 
therefore,  enables  the  identity  of  an  interacting  particle  to  be  determined,  regardless  of  its 
energy. 


Figure  3.  Time  dependence  of  scintillation  pulses  to  radiation  of  different  types  [17:227] 
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Spectrum  Unfolding 

The  organic  scintillation  detector  converts  incident  neutrons  and  gammas  to  recoil 
protons  and  Compton  recoil  electrons,  respectively.  It  is  the  energy  spectra  of  these  recoil 
protons  and  electrons  which  are  then  recorded  in  the  laboratory  by  the  detection  system. 
The  process  which  attempts  to  reproduce  the  incident  neutron  and  gamma  spectra  from 
the  measured  recoil  spectra  is  called  unfolding.  To  understand  the  unfolding  process,  one 
must  first  understand  the  relationship  between  the  incident  neutron  and  gamma  energy  and 
the  resulting  recoil  spectra.  After  that,  one  can  proceed  to  discuss  the  mathematical 
process  by  which  the  incident  spectrum  is  extracted  from  the  recoil  spectrum. 

As  stated  previously,  neutrons  entering  the  scintillator  interact  by  scattering  with 
hydrogen  and  carbon  nuclei.  The  amount  of  energy  a  neutron  can  lose  in  a  single 
collision,  and  thus  impart  to  the  recoil  nucleus,  is  dependent  upon  the  mass  of  the  target 
nucleus.  The  recoil  energy  £/{,is  given  by  [  1 7: 53 1  ] 

£,=7-^(a>s’»K  (I) 

(1  +  ^) 

where  0  is  the  scattering  angle  in  the  laboratory  reference  frame,  A  is  the  mass  number  of 
the  target  nucleus,  and  E„  is  the  energy  of  the  incident  neutron.  If  the  target  nucleus  is 
initially  at  rest  at  the  origin  of  the  lab  frame  coordinates  and  the  incident  neutron  is 
traveling  from  left  to  right  along  the  negative  x-axis,  the  scattering  angle  0  is  measured 
between  the  positive  x-axis  and  the  outbound  direction  of  the  recoil  nucleus.  To  maximize 
the  energy  deposited  in  the  course  of  the  fewest  interactions,  a  nucleus  with  a  mass  close 
to  that  of  the  neutron  is  the  preferred  target.  Since  the  mass  of  a  hydrogen  nucleus  and  a 
neutron  differ  by  less  than  0.2%,  hydrogen  makes  an  ideal  target.  It  is  possible  for  the 


10 


neutron  to  give  up  all  of  its  energy  to  the  recoil  nucleus  in  a  single  collision  with  hydrogen 
when  0=0°,  denoted  by  This  happens  when  the  incident  neutron  stops  and  the 

recoil  nucleus  continues  along  the  original  path  of  the  incident  neutron.  In  the  case  where 


0=0°  recoil,  Eq  (1)  simplifies  to  [17  531]: 


4A 


The  minimum  recoil  energy,  on  the  other  hand,  occurs  when  the  0  =90°  and  no 
energy  is  imparted  to  the  target.  With  hydrogen  as  a  target  and  assuming  isotropic 
scattering  in  the  center-of-mass  reference  frame,  the  energy  spectrum  of  the  recoil  proton 
is  ideally  rectangular,  as  shown  in  Figure  4.  In  this  figure,  the  probability  of  a  neutron 
creating  a  recoil  proton  with  a  given  energy,  P(Ep),  is  plotted  vs.  recoil  proton  energy,  Ep. 
The  maximum  energy  that  n  recoil  proton  can  receive  is  equal  to  the  energy  of  the  incident 
neutron,  E„.  This  graph  is  a  probability  distribution  function,  so  the  area  under  the  curve 


The  probability  of  a  neutron  creating  a  recoil  proton  of  a  given  energy,  P(E^,  can  also  be 
expressed  as  the  change  in  the  number  of  recoil  protons  with  changing  recoil  proton 
energy,  dN/dE. 

In  a  real  organic  scintillator,  however,  the  ideal  rectangular  response  of  Figure  4  is 
complicated  by  the  following  factors: 

1.  scattering  with  carbon, 

2.  finite  detector  resolution, 

3.  nonlinear  light  output  with  energy, 

4.  edge  effect, 

5.  multiple  scattering  with  hydrogen,  and 

6.  competing  reactions  at  high  neutron  energies  [17:538-542]. 

Carbon,  the  other  constituent  of  organic  scintillators  (xylene  is  45  atom%  carbon), 
complicates  things  somewhat.  Carbon  has  a  mass  number  of  12,  so  by  Eq.  (2)  a  neutron 
can  deposit  up  to  28%  of  its  kinetic  energy  in  a  single  elastic  scatter  with  carbon.  Thus, 
the  maximum  energy  of  a  subsequent  recoil  proton  from  hydrogen  after  a  collision  with 
carbon  varies  between  72  and  100%  of  the  original  neutron  energy.  The  energy  deposited 
during  the  carbon  scatter  is  not  converted  to  light  in  the  scintillator,  so  the  incident 
neutron  appears  to  have  a  lower  kinetic  energy  than  it  actually  had  upon  entering  the 
scintillator.  The  effect  of  carbon  scattering  is  shown  in  Figure  5.  [1 7:540] 

The  finite  resolution  of  the  detector,  due  to  sources  of  dispersion  such  as 
nonuniform  light  collection  and  photoelectron  statistics,  tends  to  wash  out  some  of  the 
distinct  structure  expected  in  the  response  function.  This  distorts  the  sharp  edge  of  the 
distribution  function  at  its  maximum  energy,  Ep,  as  shown  in  Figure  6  [17:541]. 


Figure  5.  Distortion  to  the  recoil  proton  response  function  due  to  scattering 
with  carbon.  The  thin  line  represents  the  ideal  rectangular  distribution,  while 
the  thick  line  shows  the  effect  of  carbon  scattering  [17:540] 


Figure  6,  Distortion  to  the  recoil  proton  response  function  due  to  finite  detector 
resolution.  The  thin  line  represents  the  ideal  rectangular  distribution,  while 
the  thick  line  shows  the  effect  of  finite  resolution  [17:540] 


The  light  output,  H,  from  most  organic  scintillators  is  not  linear  with  neutron 
energy,  but  is  proportional  to  where  E  is  the  energy  of  the  incident  neutron  and  k  is 


a  constant  of  proportionality.  This  relation  between  H  and  E  is  shown  in  Figure  7(a).  The 
change  in  light  output  with  changing  energy,  then,  is  the  derivative  of  this  relation 
[17:538]: 


dE 


(3) 


which  is  shown  in  Figure  7(b).  The  distortion  in  the  ideal  recoil  proton  spectrum  due  to 
nonlinear  response  is  expressed  as  the  change  in  number  of  recoil  protons  with  pulse 
height,  which  is  [17:538]: 

dN  _  dN/dE  _  constant  . 

dH/dE  ~ 

2 


Graphically,  Eq.  (4)  represents  the  division  of  Figure  4  by  Figure  7(b),  which  results  in 
Figure  7(c). 


Figure  7.  (a)Pulse  height  vs  neutron  energy;  (b)Derivative  of  (a);  (c)Distortion  of 
rectangular  proton  recoil  spectrum  due  to  nonlinear  response  of  scintillator  [17:539] 
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Organic  scintillators,  however,  do  have  a  linear  light  output  with  increasing  energy  for 
gamma  events.  For  this  reason,  spectrometry  systems  such  as  the  one  constructed  for  this 
research  are  energy  calibrated  with  a  gamma  source  of  known  energy,  usually  22Na 
[17:222], 

Another  factor  which  distorts  the  ideal  rectangular  recoil  spectrum  in  an  organic 
scintillator  is  edge  effect.  Edge  effect  occurs  when  either  a  small  detector  or  high  energy 
incident  radiation  results  in  incident  particles  exiting  the  detector  after  having  deposited 
only  part  of  their  energy  in  the  sensitive  volume  of  the  scintillator.  This  causes  a  shift  of 
highly-energetic  events  toward  lower  pulse  heights.  The  effect  on  the  response  function  is 
to  further  increase  the  slope  of  the  spectrum  in  Figure  7(c)  [17:539]. 

An  additional  effect  which  distorts  the  response  function  is  multiple  scattering 
from  hydrogen.  In  the  ideal  case,  all  of  the  energy  from  the  neutron  is  transferred  to  a 
recoil  proton  in  a  single  scatter.  However,  in  practice  it  is  possible  for  a  neutron  to  scatter 
more  than  once  with  hydrogen  nuclei  before  coming  to  rest  in  the  detector  or  escaping 
from  the  scintillator.  Because  all  such  events  usually  happen  within  a  very  short  time 
compared  to  the  fluorescence  time  of  the  scintillator,  the  light  from  multiple  scatters  is 
summed  into  a  single  pulse  whose  amplitude  is  proportional  to  the  total  light  output.  The 
effect  on  the  response  function  is  to  add  events  at  large  amplitudes  at  the  expense  of  those 
at  lower  pulse  heights,  as  shown  in  Figure  8  [17:539].  This  figure  represents  the  output 
from  a  multi-channel  analyzer  (MCA),  which  is  used  in  the  laboratory  to  view  the  pulse 
height  spectrum  of  recoil  protons.  In  this  figure,  dN/dH  is  equivalent  to  the  number  of 
counts  per  MCA  channel,  while  the  pulse  height  corresponding  to  the  neutron,  H,  is 
equivalent  to  the  channel  number  of  the  MCA. 
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Figure  8.  Distortion  tr  the  recoil  proton  response  function  due  to 
multiple  scatters  with  hydrogen  (dotted  line)  [17:541] 


As  a  final  factor,  once  the  neutron  energy  exceeds  about  8  MeV,  two  additional 
competing  reactions  must  be  considered  in  the  overall  response  of  the  scintillator.  These 
reactions  are;  [17:542] 


and 


n+  12c  +  6.17MeV->a  +  9Be 


n+  <2C  +  7.98  MeV^n'  +  3a 


The  effect  of  these  alpha-producing  reactions  is  to  increase  the  detection  efficiency  of 
neutrons  above  the  energy  thresholds  of  these  reactions  [17:542], 

The  discussion  up  to  this  point  has  not  coi.jidered  the  effect  on  the  recoil  proton 
spectrum  due  to  many  "eutrons  of  different  energies.  As  discussed  previously,  in  the  ideal 
case  a  recoil  proton  has  an  equal  probability  of  receiving  anywhere  fi’om  none  to  all  of  a 
scattered  neutron's  kinetic  energy.  Ideally  (ignoring  the  six  distortions  just  discussed),  the 
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measured  recoil  proton  spectrum,  M(H^,  resulting  from  several  monoenergetic  neutrons 
will  match  the  probability  distribution  function  of  Figure  4.  For  the  same  number  of 
neutrons  of  any  different  neutron  energy,  we  expect  a  rectangular  response  with  a  different 
height,  but  with  the  same  area  under  the  curve.  If  we  were  to  measure  neutrons  of  several 
discrete  energies,  each  discrete  energy  group  having  the  same  number  of  neutrons,  we 
would  ideally  expect  a  family  of  rectangular  responses,  R{HpJE„),  each  with  the  same  area, 
as  shown  in  Figure  9. 


Fi^re  9.  Ideal  response  curves  of  neutrons  of  several  discrete  energies 


Therefore  the  response,  R{Hp,E„),  is  the  probability  that  a  neutron  with  energy  E„  will 
scatter  once  with  a  hydrogen  nucleus  and  create  a  recoil  proton  of  energy  Ep,  which  will  in 
turn  produce  a  pulse  of  height  Hp  in  the  detection  system.  If  every  neutron  that  enters  the 
scintillator  interacts,  then  the  following  relationship  holds  for  each  response  curve  [5:85]; 

=  l  (5) 
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The  response  function,  R{Hp^„),  is  ideally  l/E„  for  0<Ep^En  [5;85]. 

If  the  curves  from  Figure  9  are  added  together  to  simulate  incident  neutrons  of 
several  discrete  energies,  a  plot  of  the  total  number  of  recoil  protons,  for  all  of  the 

discrete-energy  neutrons  that  entered  the  scintillator  volume  is  obtained,  as  shown  in 
Figure  10. 


The  measured  recoil  proton  spectrum  is  related  to  the  true  incident  neutron 

spectrum,  as  follows  [17:672]: 

^  =  /'■"  liH,.  E,)nE,)dE.  (6) 

where  R{Hp^„)  is  called  the  response  function  or  response  matrix  and  is  specific  to  the 
spectrometry  system  used.  The  response  matrix  accounts  for  the  response  of  the  detector 
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to  incident  neutrons  of  each  energy  E„  in  terms  of  the  pulse  height  of  the  corresponding 
recoil  proton,  Hp.  Taken  together,  the  term  in  the  integral  is  the  differential  probability 
that  a  quantum  of  energy  within  dE„  about  E„  leads  to  a  pulse  with  amplitude  within  dH 
about// [17:672]. 

The  summed  response  functions  shown  in  Figure  10  still  assume  that  the  response 
of  the  detector  to  neutrons  of  each  discrete  energy  is  ideal  and,  therefore,  rectangular.  If 
all  of  the  distortions  to  the  ideal  rectangular  response  are  now  considered,  the  recoil 
proton  spectrum  of  a  multi-energy  neutron  source,  such  as  PuBe,  becomes  quite 
complicated.  The  response  function  of  the  detector  used  must  be  known  with  a  high 
degree  of  accuracy  to  unfold  recoil  proton  spectra  and  correctly  reproduce  the  incident 
neutron  spectrum  from  the  measured  recoil  proton  spectrum.  For  this  reason,  the 
response  functions  of  a  particular  detector  to  monoenergetic  neutrons  of  several  different 
energies  are  either  carefully  measured  in  the  laboratory  or  are  calculated  using  a  Monte 
Carlo  code  in  which  the  detector  has  been  carefully  modeled.  Verbinski  et  al.  [29:8] 
measured  the  response  of  an  NE-213  spectrometry  system  to  monoenergetic  neutrons  at 
20  energies  between  0.2  aiid  20  MeV  and  compared  these  measurements  to  Monte  Carlo 
calculations.  These  results  form  the  basis  for  the  neutron  response  matrix  used  by  the 
FORIST  unfolding  code.  For  this  reason,  many  neutron  spectrometry  experiments  are 
now  performed  using  the  laboratory  protocols  of  Verbinski  et  al.  in  order  that  FORIST 
and  its  neutron  response  matrix  may  be  utilized. 

The  integral  unfolding  problem  given  in  Eq  (6)  must  be  reduced  to  a  matrix 
unfolding  problem  to  solve  it  numerically.  The  matrix  equation  is  [5:85;17:673]: 

w(h,),=Z«(w,, £.)/{£.),  (7) 

J 
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where  A/,  is  the  measured  number  of  counts  in  the  channel  of  a  multi-channel  a 
Rij  is  the  response  matrix  coupling  the  pulse  height  interval  with  the energy 
and  Tj  is  the  true  neutron  intensity  in  the energy  interval  [  1 7:673],  M and  Tar 
vectors,  while  7?  is  a  square  matrix  [5:85],  The  matrix  R  is  upper  triangular  becai 
recoil  proton  can  not  receive  more  energy  than  the  incident  neutron  had  before  th 
collision. 

Eq  (7)  can  be  solved  for  the  incident  neutron  spectrum,  T,  by  multiplying 
sides  of  the  equation  by  the  inverse  of  7?  [15:4]: 

i 

Problems  occur,  however,  in  the  numerical  solution  because  of  round-off  tti  c  s . 
the  matrix  inversion  process,  as  well  as  uncertainties  in  the  detector  readings  and  i;i  the 
response  matrix.  These  errors  can  produce  large  oscillations  about  the  true  spectrum  in 
the  solution  vector.  In  practice,  errors  in  R  and  M  are  inevitable,  and  the  more  channels  of 
data  that  are  used,  the  more  likely  it  is  for  R  is  to  be  ill-conditioned  and  the  system  to 
become  unstable,  so  that  small  errors  produce  large  changes  in  the  output  spectrum 
[8:131].  Most  unfolding  methods,  therefore,  must  have  some  means  of  limiting  errors  and 
oscillations.  Some  strategies  include  forcing  the  solution  vector,  T,  to  be  "near  z  ' 

[15:5],  or  to  have  a  certain  degree  of  smoothness  [8:131], 

The  unfolding  of  the  incident  gamma  ray  spectrum  from  the  measured  C(  :  . 
recoil  electron  spectrum  proceeds  in  a  very  similar  manner.  The  interaction  pror ; . 
Compton  scattering  takes  place  between  the  incident  gamma  photon  and  an  elec  .  .  '  i 

scintillator.  In  Compton  scattering,  the  incoming  gamma  is  deflected  through  ai  .  4 
(in  the  lab  frame)  with  respect  to  its  original  direction.  The  gamna  transfers  a  p  . .  '  : ;  :  t' 
its  energy  to  the  electron,  which  is  assumed  to  be  at  rest  initially.  The  energy  tn  ii'  rred 
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from  the  gamma  to  the  recoil  electron  is  dependent  upon  the  scattering  angle  as  follows 


[l-;53]; 


— ^(l-cos^) 
m^c 

£  ^ 

1  +  —-^(1-  cos 


where  Ey  is  the  energy  of  the  incident  gamma,  Eg  is  the  energy  of  the  resulting  recoil 
electron,  and  nigC^  is  the  rest-mass  energy  of  the  electron  (0.5 1 1  MeV).  The  maximum 
energy  transferred  to  the  electron  occurs  when  0=180°.  Under  normal  circumstances,  all 
scattering  angles  will  occur  in  the  detector.  Therefore  a  continuum  of  energies  can  be 
transferred  to  the  electron,  ranging  from  zero  up  to  the  maximum  predicted  by  Eq  (9) 
[17:290-291].  For  a  monoenergetic  gamma  say,  the  recoil  electron  energy  distribution 
with  scattering  angle  has  the  shape  shown  in  Figure  1 1 . 

Note  that  there  is  some  similarity  between  the  energy  distribution  of  recoil 
electrons  in  Figure  1 1  and  the  idea!  distribution  of  recoil  protons  in  Figure  4.  Similarly, 
too,  the  ideal  energy  distribution  of  the  recoil  electrons  suffers  distortions  due  to  factors 
such  as  the  finite  resolution  of  the  detector.  One  important  feature  that  the  recoil  electron 
spectrum  does  not  experience  is  the  nonlinear  light  collection  with  increasing  energy  | 

1 

observed  with  the  recoil  proton  spectrum.  This,  coupled  with  the  fact  that  the  Compton 
edges,  though  smeared  by  the  finite  resolution  of  the  detector,  are  identifiable  in  the  recoil 
electron  spectrum,  makes  gamma  ray  sources  ideal  for  energy  calibrating  organic  1 
scintillation  detectors.  The  exact  theory  and  process  of  energy  calibration  are  described  in 
following  sections. 
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Figure  11.  Energy  distribution  of  Compton  recoil  electrons  [17:291] 


As  with  neutrons,  the  response  of  a  detector  to  incident  gamma  rays  can  be  either 
measured  or  calculated,  yielding  a  gamma  response  matrix  for  the  detector  which  allows 
the  recoil  electron  spectra  to  be  unfolded  to  reproduce  the  incident  gamma  spectra. 
Ingersoll  and  Wehring  [13:551]  used  the  NE-213  detector  configuration  of  Verbinski  et 
al.  to  measure  the  gamma  ray  response  functions  of  the  spectrometry  system  for  100 
monoenergetic  gamma  rays  from  2  to  11.5  MeV.  These  measurements  were  used  to 
construct  the  gamma  response  matrix  which  is  used  with  the  FORJST  code. 


Energy  Calibration 

To  unfold  a  recoil  spectrum  to  yield  the  incident  spectrum,  the  spectrometer  must 
be  energy-calibrated.  Calibration  is  required  both  to  assign  an  energy  scale  to  the 
unfolded  spectrum  and  to  match  the  detector's  response  to  the  response  matrices  provided 
with  FORIST. 


^saagg 


The  pulse  height  response  of  the  scintillator  is  nonlinear  for  incident  neutrons,  but 
it  is  linear  for  incident  gammas,  so  the  easiest  way  to  calibrate  the  detection  system  is  to 
use  a  gamma  source  and  to  perform  a  linear  fit  to  determine  the  calibration.  A  gamma 
source  with  at  least  two  peaks  of  known  energy  is  required.  Following  the  protocols  of 
Verbinski  etal.  [29:1*1],  a  ^^Na  source,  which  has  two  gammas  of  energy  0.51 1  and  1.275 
MeV,  is  used  to  calibrate  the  system.  The  maximum  Compton  recoil  electron  energies 
produced  by  these  two  gammas,  calculated  using  Eq  (9),  are  0.341  and  1.062  MeV 
respectively.  A  22Na  Compton  recoil  electron  spectrum  from  the  literature  is  shown  in 
Figure  12. 


Figure  12.  Compton  recoil  electron  spectrum  from  literature  [18:523] 
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Verbinski  et  al.,  in  their  calibration  work,  defined  an  energy  unit  called  a  "light 
unit"  (!u).  They  defined  0.885  lu  to  be  equal  to  the  maximum  recoil  energy  of  the  1.275- 
MeV  gamma  from  22Na,  or  equal  to  1 .062  MeV.  [6: 183]  The  FORJST  unfolding  code 
requires  that  the  recoil  spectra  recorded  in  the  lab  be  calibrated  in  terms  of  light  units. 

Because  of  the  finite  resolution  of  the  detector,  the  Compton  edge  is  not  a  sharp 
event,  but  is  smeared  over  several  MCA  channels,  as  can  be  seen  in  Figure  12.  Various 
researchers  have  differed  in  their  choice  of  the  location  of  the  maximum  recoil  electron 
energy  on  these  Compton  edges,  selecting  either  the  peak,  two-thirds  height,  or  half¬ 
height  [7:351;18;522].  The  half-height  is  used  for  this  thesis  research  because  it 
I 

consistently  provides  more  accurate  results  [24;267;16;341]. 

The  calibration  is  calculated  as  follows.  The  channel  locations  of  the  half-heights 
of  the  two  Compton  edges  are  used  to  define  a  line.  First,  the  slope  of  the  line  is 
calculated  using  the  equation: 

w  =  (10) 

Xj-JC, 

where  m  is  the  slope  of  the  line,  (xj  j'j)  are  the  half-height  channel  location  and  maximum 
recoil  electron  energy  of  the  first  (0.5 1 1  MeV)  gamma,  and  (X2^y2)  are  the  half-height 
channel  location  and  maximum  recoil  electron  energy  of  the  second  gamma  (1 .275  MeV). 
Next,  the  y-intercept  of  the  line  is  calculated  using  the  equation: 

* 

b  =  y^-mx^  (11) 


where  h  is  the  y-intercept  of  the  line  calculated  in  Eq  (1 1),  is  the  slope  calculated  in 


Eq  (10),  aiid  (x2j'2)  again  correspond  to  the  half-height  channel  location  and  maximum 
recoil  electron  energy  of  the  second  gamma.  Third,  the  zero-intercept,  which  is  the 
location  where  zero  energy  crosses  the  x-axis,  is  found  using  the  equation: 


(12) 


where  z  is  the  zero-intercept,  while  m  and  h  are  the  slope  and  y-intercept  calculated 
previously.  Finally,  the  gain,  in  light  units,  is  calculated  for  the  system  using  the  equation: 


G  = 


0.885 


(13) 


where  G  is  the  gain,  z  is  the  zero  intercept  from  Eq  (12),  Xi  is  the  channel  location  of  the 
maximum  recoil  energy  of  the  second  gamma,  and  0.885  is  the  light  unit  conversion  factor 
of  Verbinski  et  al.  Detailed  discussions  of  the  calibration  procedure  and  the  gain  and 
zero-intercept  requirements  of  the  FORIST  code  are  presented  the  "Laboratory 
Procedure"  section. 
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Ill.  Description  of  Apparatus 


The  detector  used  in  this  research  project  is  a  xylene-based  organic  scintillator 
similar  to  the  one  used  by  Verbinski  et  al.  This  type  and  size  of  detector  were  chosen  so 
that  the  FORI  ST  unfolding  code  could  be  used  to  reproduce  the  incident  neutron  and 
gamma  spectra  from  the  recorded  recoil  proton  and  Compton  electron  spectra.  This 
section  first  compares  the  detector  used  in  this  research  to  the  one  used  by  Verbinski  ^ 
et  al.  to  show  that  the  FORIST  response  matrices  apply  to  the  current  detector.  Next,  it 
presents  the  electronics  used  in  the  laboratory  to  record  the  data,  with  an  emphasis  on  the 
setup  used  to  differentiate  between  neutron  and  gamma  events  on  the  basis  of  PSD. 
Following  the  discussion  of  the  electronics,  the  requirement  for  and  design  of  the  shadow 
shield  are  discussed.  Next,  the  radiation  sources  used  during  the  course  of  this  research 
are  described.  Finally,  this  section  presents  a  brief  description  of  all  of  the  computer 
software  used  in  conjunction  with  the  spectrometry  system,  including  the  FORIST 
unfolding  code. 


Detector 

The  detector  used  in  this  research  consists  of  a  5.0-cm  by  5.0-cm  cylindrical 
Bicron(BC)-501  A  liquid  organic  scintillator  coupled  to  a  Hamamatsu  R329-02 
photomultiplier  tube  (PMT)  by  means  of  a  1.6-cm  thick  clear  plastic  light  pipe.  Both  the 
scintillator-light-pipe  interface  and  the  light-pipe-PMT  interface  are  optically  coupled  with 
a  Dow  Coming  transparent  silicon  fluid  which  has  a  kinematic  viscosity  of 2000  stokes. 
Aluminum  foil  is  used  as  a  reflective  coating  on  the  outside  of  the  light  pipe,  and  the  entire 
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assembly  is  sealed  with  black  electrical  tape  to  eliminate  light  leaks.  The  PMT  was  then 
surrounded  by  a  p-metal  shield  to  reduce  magnetic  interference.  An  equipment  list  is 
included  in  this  thesis  in  Appendix  A,  while  a  diagram  of  the  detector  is  shown  in 
Figure  13. 


Figure  13.  Diagram  of  detector  used  in  this  research  project 


The  detector  used  by  Verbinski  et  al.  consisted  of  a  4.60-cm  by  4,65-cm  diameter 
NE-213  liquid  organic  scintillator  mounted  on  a  0.63-cm  thick  clear  plastic  light  pipe  and 
coupled  to  an  RCA-6810-A  PMT.  An  aluminum  foil  reflector  was  used  to  cover  the 
scintillator  and  light  pipe  [6;  183]. 


The  NE-213  liquid  scintillator  manufactured  by  Nuclear  Enterprises  and  used  in 
neutron  spectrometry  for  many  years  is  no  longer  available  in  the  United  States.  An 
available  substitute,  and  the  one  used  in  the  detection  system  constructed  for  this  project, 
is  the  BC-501 A  liquid  scintillator  manufactured  by  Bicron.  A  comparison  of  these  two 
scintillators,  shown  in  Table  1,  shows  that  they  are  essentidly  identical. 


Table  1.  Comparison  ofNri-213  and  BC-501A  scintillators  [17:221;3:11 


Property 

NE-213 

BC-501A 

Light  output,  %  anthracene 

78 

78 

Density  (g/cm^) 

0.874 

0.874 

Refractive  Index 

1.508 

1.505 

Wavelength  of  max  emission  (run) 

425 

425 

H/C  (atomic  ratio) 

1.213 

1.212 

The  cylindrical  geometry  of  the  NE-213  and  BC-501  A  scintillators  are  the  same, 
but  there  may  be  a  small  difference  in  the  dimensions  of  the  cylinder.  The  5.0  by  5.0-cm 
dimensions  reported  for  the  scintillator  used  in  the  detection  system  in  this  thesis  are  the 
outer  dimensions  of  the  plastic-and-aluminum  container  in  which  the  liquid  scintillator  is 
contained.  Verbinski  et  al.  do  not  mention  whether  their  reported  dimensions  are  the 
inner  or  outer  measurements  of  their  scintillator,  but  judging  from  the  number  of  other 


researchers  that  have  reported  using  a  5.0  by  5.0-cm  scintillator  in  conjunction  with  the 
FORIST  code  [13:55 1;16;337],  either  Verbinski  e/L/.  probably  reported  the  inner 
dimensions  of  their  scintillator. 

The  light  pipe  used  in  the  current  detector  is  the  same  as  Verbinski's  except  that  it 
is  0.97  cm  thicker  than  his.  This  slightly  thicker  dimension  was  chosen  to  allow  the  rim  of 
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the  BC-501 A  scintillator  to  sit  directly  on  top  of  the  |i-meta!  shield  without  cutting  0.97 
cm  from  the  top  of  the  shield.  The  light  pipe  serves  only  to  distribute  the  scintillator  light 
more  uniformly  on  the  photocathode  of  the  PMT ,  so  the  absence  of  a  light  pipe  might  be 
expected  to  affect  the  response,  but  the  small  difference  in  thickness  should  have  a 
negligible  effect  on  the  response  of  the  detector  to  incident  radiation  [1 7:242-246]. 

Finally,  one  must  compare  the  parameters  of  the  two  photomultiplier  tubes  used, 
the  RCA  68 1 0-i'j.  of  Verbinski  et  al.  and  the  Hamamatsu  R329-02  of  the  current  detection 
system.  The  Hamamatsu  tube  used  with  the  current  detector  configuration  was  marketed 
by  Bicron  as  the  properly-matched  PMT  for  the  BC-501  A  scintillator.  The  RCA  86 10- A 
tube  used  by  Verbinski  et  al.  is  no  longer  made,  but  the  tube  marketed  by  RCA  as  its 
replacement,  the  RCA  8575,  has  similar  characteristics.  The  RCA  and  the  Hamamatsu 
tubes  are  compared  in  Table  2.  Both  Verbinski's  PMT  and  the  tube  used  for  this  project 
are  very  well  matched  to  the  425-nm  wavelength  of  maximum  emission  of  the  scintillator, 
and  both  have  response  times  appropriate  to  PSD.  The  selection  of  the  Hamamatsu  tube 
should  allow  the  use  of  Verbinski's  response  functions  with  this  detector. 

Table  2.  Comparison  of  RCA  8610-A  and  Hamamatsu  R329-02 


ihotomultiplier  tubes  [27:12-13;!  1:44-45 


Parameter 

RCA  8610-A  (8575) 

Hamamatsu  R329-02 

Spectral  Range  (nm) 

260  -  650 

300  -  650 

Peak  Wavelength  (run) 

420 

420 

Voltage  Divider  Stages 

12 

12 

Maximum  HV  (V) 

3000 

2700 

Diameter  (mm) 

51 

51 

Rise  Time  (ns) 

2.8 

2.6 

Transit  Time  (ns) 

37 

48 
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Electronics 


The  electronics  used  in  the  spectrometry  system  constructed  for  this  thesis  can  be 
divided  into  four  subsystems;  the  detection  subsystem,  the  PSD  subsystem,  the  energy 
subsystem,  and  the  output  subsystem.  The  function  of  each  of  the  four  subsystems  is 
discussed  in  turn  though  a  detailed  description  of  the  operation  of  each  electronic  module 
is  described  in  tne  "Experimental  Procedure"  section.  A  diagram  of  the  spectrometry 
system's  electronics  is  presented  in  Figure  14,  and  an  equipment  list  is  provided  in 
Appendix  A. 


Figure  14.  Block  diagram  of  electronic  equipment  for  spectrometry  system 


The  function  of  the  electronics  in  the  detection  subsystem  is  to  provide  the 
optimum  signal  from  the  scintillator  for  the  PSD  subsystem  and  the  energy  subsystem. 

This  requires  a  fast,  unsaturcted  signal  to  convey  the  PSD  information  and  to  properly 
convey  the  energy  information.  The  electronics  in  the  detection  subsystem  consist  of  the 
photomultiplier  tube  base,  high  voltage  (HV)  supply,  and  preamplifier. 

The  PMT  base  accepts  the  input  from  the  HV  supply  and  has  BNC  outputs  for 
bodi  the  anode  and  the  9**’  dynode  signal  from  the  PMT.  Many  early  researchers, 
including  Verbinski  el  ai,  took  the  output  from  a  dynode  for  their  energy  signal  and  the 
output  from  both  a  later  dynode  and  the  anode  for  their  PSD  signal.  This  method  usually 
required  modifications  to  the  PMT  base  and  a  number  of  complicated  modules  in  the  PSD 
subsystem  [29: 13].  This  is  no  longer  necessary  since  nuclear  electronics  have  progressed 
since  the  original  work,  however.  The  electronic  configuration  chosen  for  the  current 
research  project  was  based  on  a  setup  suggested  in  the  operating  manual  for  an  ORTEC 
model-460  delay-line  amplifier  [10: 10].  In  this  new  setup,  both  the  timing  and  energy 
signals  are  taken  from  the  9^  dynode.  The  advantages  of  this  setup  are  that  the  PMT  base 
does  not  have  to  be  modified  and  the  output  is  suitable  for  processing  by  simple  sets  of 
follow-on  electronics  in  the  timing  and  energy  subsystems. 

The  HV  supply  provides  the  potential  difference  to  the  PMT  necessary  to 
accelerate  the  electrons  from  one  dynode  to  the  next.  The  HV  supply  must  be  chosen  to 
properly  match  both  the  magnitude  and  polarity  of  the  PMTs  voltage  requirements.  For 
this  system  the  HV  supply  provides  negative  high  voltage  to  the  photocathode. 

The  preamplifier  serves  as  an  intermediate  amplification  step  in  the  system,  which 
improves  the  signal-to-noise  ratio  and  minimizes  capacitive  loading  on  the  detector.  The 
connection  from  the  9**’  dynode  of  the  PNff  base  to  the  preamp  must  be  transmitted  on 
short  a  cable  to  minimize  degradation  of  the  signal.  The  output  from  the  preamp  can  be 
transmitted  over  long  lengths  of  930  coaxial  cable  to  allow  the  operator  to  locate  the  rest 
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of  the  electronics  at  a  safe  distance  from  the  radiation  sources.  During  this  research,  the 
operator's  station  was  located  approximately  4  meters  from  the  detector  and  radiation 
sources. 

The  next  subsystem  to  be  discussed  is  the  PSD  subsystem.  The  function  of  the 
pulse  shape  discrimination,  or  timing,  subsystem  is  to  process  the  available  pulse  shape 
information  produced  by  the  organic  scintillator  and  to  use  this  information  to  establish  a 
gate  to  select  gamma  or  neutron  events.  The  PSD  subsystem  chosen  for  this  research 
consists  of  a  delay  line  amplifier  (DLA),  a  pulse  shape  analyzer  (PSA),  and  a  gate  and 
delay  generator  (GDG). 

The  delay  line  amplifier  takes  the  positive  unipolar  signal  from  the  preamp  and 
produces  a  double-delay-line  shaped  bipolar  pulse  with  positive  polarity  leading  [10;  1]. 
Timing  measurements  are  obtained  from  the  zero  crossover  point  of  the  bipolar  pulse. 
Neutrons,  which  have  a  larger  delayed  fluorescence  component  in  the  scintillator  than 
gammas,  will  produce  pulses  in  the  DLA  which  cross  zero  later  than  those  pulses 
produced  by  gamma  events.  A  sample  output  from  the  DLA  is  shown  in  Figure  20  of  the 
"Experimental  Procedure"  section. 

The  pulse  shape  analyzer,  the  next  module  in  the  PSD  subsystem,  measures  the 
fall  time  of  an  input  pulse  and  generates  a  linear  output  pulse  with  an  amplitude  that  is 
proportional  to  the  fall  time.  A  logic  output  pulse  is  also  generated  when  the  fall  time  of 
the  input  pulse  lies  between  the  upper-  and  lower-level  discriminators  of  the  PSA  time 
window  [9;  1  ].  In  the  context  of  this  spectrometry  system,  the  PSA  operates  on  the 
positive  lobe  of  the  bipolar  output  pulse  from  the  DLA.  The  linear  output  pulse  generated 
is  suitable  for  input  to  an  analog-to-digital  converter  (ADC)  and  multi-channel  analyzer 
(MCA).  This  linear  pulse  provides  a  timing  spectrum  on  the  MCA  which  can  be  used  to 
set  the  PSA  time  window  to  trigger  a  logic  pulse  only  when  a  gamma  event  is  detected  by 
the  scintillator.  This  logic  pulse  is  then  used  as  the  basis  for  the  gate  between  neutron  and 

32 


gamma  events.  Sample  outputs  from  the  PSA  are  shown  in  Figures  21  and  22  of  the 
"Experimental  Procedure"  section. 

The  gate  and  delay  generator,  the  last  module  in  the  PSD  subsystem,  is  used  to 
shape  the  logic  pulse  output  from  the  PSA  in  both  amplitude  and  width  to  meet  the 
specifications  required  by  the  ADC  for  a  gate  pulse.  The  logic  output  from  the  PSA  is 
used  as  the  input  to  the  GDG,  and  the  output  from  the  GDG  is  used  as  the  gate  input  to 
the  ADC,  which  is  the  first  module  in  the  output  subsystem.  Sample  outputs  from  the 
GDG  are  shown  in  Figures  23  and  24  of  the  "Experimental  Procedure"  section. 

The  third  subsystem  to  be  described  is  the  energy  subsystem.  The  function  of  the 
energy  subsystem  is  to  take  the  energy  information  from  the  detector  and  process  it  to 
meet  the  specifications  required  by  the  ADC  for  a  signal  pulse.  The  energy  subsystem 
consists  of  an  amplifier  and  a  delay  amplifier. 

The  amplifier  increases  the  amplitudes  of  the  pulses  from  the  preamp  to  a 
magnitude  suitable  for  the  ADC.  The  gain  of  the  amplifier  is  selected  to  match  the  gain 
settings  required  for  the  FORIST  unfolding  code.  The  specific  gain  settings  are  described 
later.  A  sample  output  from  the  amplifier  is  shown  in  Figure  24. 

The  delay  amplifier  is  used  to  match  a  pulse  from  the  energy  subsystem  with  the 
corresponding  pulse  from  the  PSD  subsystem.  Because  there  are  more  modules  in  the 
PSD  ^bsystem,  the  timing  signal,  without  some  way  to  delay  the  energy  pulse,  would 
reach  the  ADC  later  than  the  energy  pulse.  A  delay  of  a  few  microseconds  is  required. 

An  osalloscope  is  used  to  view  the  outputs  from  the  GDG  and  amplifier  in  order  to  adjust 
this  delay.  The  output  of  the  delay  amp  is  then  used  as  the  signal  input  to  the  ADC.  The 
gate  and  signal  pulses  are  shown  in  Figure  24. 

j'he  fourth  and  final  subsystem  in  the  spectrometer  is  the  output  subsystem.  The 
function  of  the  output  subsystem,  which  consists  of  the  analog-to-digital  converter  (ADC) 
and  multi-channel  analyzer  (MCA),  is  to  use  the  PSD  information  to  select  only  gamma  or 
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neutron  events  and  to  then  collect  the  corresponding  Compton  electron  or  recoil  proton 
spectrum  in  the  MCA. 

The  ADC  receives  the  energy  signal  from  the  energy  subsystem  and  is  gated  by  the 
timing  gate  from  the  PSD  subsystem.  A  switch  on  the  ADC  allows  selection  of  gamma 
events,  neutron  events,  or  a  coinbtnation  of  the  two. 

The  MCA,  which  is  a  program  running  on  a  personal  computer  in  thi  :  .  )f  this 
spectrometry  system,  assigns  the  incoming  energy  pulses  to  bins,  or  channels, .  ' ;  jct 
proportion  to  the  ainplitude  of  the  signal.  The  resulting  recoil  proton  or  Conri  ■  n 
electron  spectra  recorded  by  the  MCA  are  then  saved  to  a  floppy  disk  and  later  i  nifolded 
by  the  FORIST  unfolding  code  to  generate  the  incident  neutron  and  gamma  spectra. 


To  measure  just  those  neutrons  or  gammas  traveling  directly  from  the  source  to 
the  detector,  the  background  contributions  to  the  spectrum  must  be  subtracted.  The 
background  spectrum  results  from  neutrons  and  gammas  from  the  source  scattering  off  the 
floor,  walls,  and  other  scattering  surfaces  in  the  experimental  area,  as  well  as  any  radiation 
emitted  by  any  other  radioactive  sources  which  may  be  in  the  area.  To  measure  the 
background,  or  room  return  spectrum  from  the  radiation  source,  a  shadow  shield  is  used 
to  attenuate  the  direct  source  spectrum.  The  spectrum  recorded  with  the  shadow  shield 
can  then  be  subtracted  from  the  spectrum  of  the  direct  source  recorded  without  the 
shadow  shield  to  reduce  the  room  return  contribution. 

An  effective  shadow  shield  must  completely  obstruct  the  line  of  sight  between  the 
source  and  the  detector,  and  must  be  composed  of  materials  which  are  effective  shields  for 
gammas  and  neutrons.  The  shield  assembled  for  this  spectrometry  system  consists  of  36 
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inches  of  paraffin  followed  by  4  inches  of  lead  in  a  rectangular  array  of  total  dimensions 
40"  X  11 "  X  6".  The  11 "  x  6"  rectangle  is  the  cross  section  along  the  axis  between  the 
source  and  detector,  which  is  large  enough  to  obstruct  the  line  of  sight. 

The  paraffin  used  in  the  shield,  which  is  contained  in  two  18-inch  long  plywood 
boxes,  was  chosen  because  it  is  a  hydrogenous  material  and,  thus,  and  excellent  neutron 
absorber.  The  36  inches  of  paraffin  between  the  source  and  detector  provide  13  mean  free 
paths  of  absorber  for  monoenergetic  14-MeV  neutrons,  and  even  more  MFPs  for 
neutrons  of  all  lower  energies.  This  means  that  the  shield  only  passes  0.0023%  of  incident 
14-MeV  neutrons  and  even  fewer  neutrons  of  lower  incident  energies. 

The  lead  used  in  the  shield  is  composed  of  two  2-inch  thick  lead  bricks,  '^he 
4  inches  of  lead  between  the  source  and  detector  provide  5.6  mean  free  paths  of  absorber 
for  10-MeV  gammas,  and  even  more  MFPs  for  gammas  of  all  lower  energies.  This  means  / 

that  the  shield  only  passes  0.38%  of  incident  10-MeV  ganunas  and  even  fewer  gammas  of 
lower  incident  energies. 

The  lead  is  placed  behind  the  paraffin,  as  viewed  from  the  source,  to  absorb  any 
gammas  made  by  (n,y)  events  within  the  paraffin  caused  by  the  absorption  of  the  high- 
energy  incident  neutrons  [23:57-58].  There  was  no  need  for  additional  paraffin  after  the 
lead,  as  (y,n)  reactions  in  the  lead  are  highly  unlikely. 


Radiation  Sources 

Four  radioactive  sources  were  used  during  the  course  of  this  research,  two 
239puBe  sources  and  two  22Na  sources.  The  PuBe  sources  provide  both  neutrons  and 
gammas,  while  the  22Na  provide  only  gammas.  The  four  sources  are: 

1.  PuBe  #T022:  78  mCi,  measured  22  September  1992; 
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2.  PuBe  #M-1 170;  4.7  Ci,  measured  22  September  1992; 

3.  22Na#217:  4.7  pCi,  measured  22  September  1992;  and 

4  22^3  #T028;  97.6  pCi,  measured  1  October  1992. 

The  weaker  of  the  two  PuBe  sources,  #T022,  was  used  during  all  of  the  setup  and 
optimization  of  the  spectrometer,  as  well  as  for  most  of  the  data  runs,  in  order  to  limit  the 
radiation  exposure.  The  smaller  of  the  two  22Na  sources,  #217,  was  used  as  the  gamma- 
only  source  throughout  the  setup  and  optimization  work,  as  well  as  for  nearly  all  of  the 
energy  calibration  data  runs.  The  large  PuBe  source,  #M-1 1 70,  was  used  for  only  two 
data  runs,  in  which  a  higher  count  rate  was  desired  to  improve  the  counting  statistics.  The 
large  22Na  source,  #T028,  was  used  only  once  in  a  data  run  to  measure  the  22Na  gamma 
spectrum  with  the  spectrometry  system.  The  #T028  source  was  also  used  to  calibrate  the 
spectrometer  for  that  run  only. 

The  239puBe  neutron/gamma  sources  were  used  for  two  reasons.  First,  they 
provide  neutron  and  gamma  spectra  which  are  both  well-characterized  and  within  the 
detection  range  of  the  spectrometry  system  (1  to  14  MeV).  Second,  the  sources  are  quite 
compact,  making  them  easy  to  handle,  store,  and  place  in  position  in  relation  to  the 
detector.  The  spectrum  generated  by  the  PuBe  source  is  quite  complicated,  so  a  detailed 
discussion  of  the  PuBe  source  is  presented  in  Appendix  B  rather  than  here  in  the  body  of 
the  thesis. 

The  22Na  gamma  sources  were  used  in  order  to  follow  the  convention  of  Verbinski 
et  al.  This  simple  source  has  gamma  pf»ks  at  two  energies:  0.5 1 1  MeV  and  1 .275  MeV. 


Computer  Software 


The  spectrometry  system  constructed  for  this  research  requires  the  use  of  four 
pieces  of  computer  software.  The  most  obvious,  of  course,  is  the  FORIST  unfolding 
code.  Additionally,  the  system  uses  a  program  called  ZSHIFT,  which  takes  the  recoil 
spectra  recorded  on  the  MCA,  along  with  the  appropriate  calibration  information,  and 
prepares  the  input  file  for  FORIST  to  read.  The  output  file  from  FORIST  is  interpreted  by 
a  third  program,  called  GRAPHIC  S.TK,  which  produces  plots  of  the  spectra  of  the 
incident  neutrons  and  gammas.  The  fourth  piece  of  software,  called  CALIBER.  TK, 
generates  the  calibration  information  required  by  FORIST  and  ZSHIFT.  FORIST  and 
ZSHIFT  are  both  written  in  FORTRAN  and  were  modified  to  run  on  a  VAX  6420 
mainframe  computer.  The  other  two  codes  are  written  using  a  software  package  called 
"TK  Solver  Plus"  [28;  1],  and  are  meant  to  run  on  any  IBM-compatible  personal  computer. 
The  four  codes  are  presented  in  their  entirety  in  Appendices  C,  D,  E,  and  F,  while 
instructions  for  using  the  codes  are  given  in  the  "Experimental  Procedure"  section  of  this 
thesis.  All  of  the  code  examples  provided  in  the  appendices  correspond  to  the  same  data 
run  (PuBe  #M-1 170  neutrons,  run  #9),  so  they  may  be  studied  as  a  set  to  see  how  the  four 
codes  work  together. 
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IV.  Experimental  Procedure 


This  section  is  intended  to  be  a  user's  guide  to  operating  the  spectrometry  system 
constructed  for  this  project,  or  even  for  the  construction  and  use  an  identical  system  from 
scratch.  This  section  first  describes  the  arrangement  of  tne  equipment  in  the  laboratory'. 
Next,  it  describes  how  to  set  up  and  optimize  the  electronics.  Then  it  describes,  step  by 
step,  how  to  collect  data  using  this  spectrometer.  Finally,  it  gives  instmctions  for 
processing  and  unfolding  the  laboratory  data  to  generate  the  spectrum  of  the  incident 

1 

radiation. 


Laboratory  Arrangement  i 

i 

The  detector,  source,  and  shadow  shield  (when  used)  are  placed  on  a  desk  top  ip 
the  configuration  shown  in  Figure  15.  The  source-detector  separation  distance  is 
42  inches,  or  just  2  inches  larger  than  the  length  of  the  shadow  shield.  This  minimum  I 

i 

distance  was  chosen  because  of  the  relative  weakness  of  the  small  PuBe  source  (#T022]), 
which  would  have  required  inordinately  long  counting  times  if  the  distance  had  been 
increased  any  further.  The  detector  assembly  is  placed  in  a  wooden  cradle  on  top  of  a 
paraffin-filled  1 8"  x  1 1"  x  6"  plywood  box  to  raise  the  detector  at  least  a  small  distance 
above  the  desktop.  The  detector  is  placed  so  that  the  radiation  from  the  source  is  incident 
on  the  curved  face  of  the  detector,  in  accordance  with  the  protocols  of  Verbinski  et  al 
[29:8].  The  source  under  study  is  placed  in  a  clip  on  a  ring  stand  located  at  the  source 
location  and  raised  to  the  same  height  above  the  desktop  as  the  detector.  Between  the 
source  and  detector  (below  the  source-detector  plane)  are  two  more 
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Figure  1 5.  Laboratory  arrangement  of  source,  detector,  and  shadow  shield 


paraffin-filled  plywood  boxes,  on  which  the  shadow  shield  rests  when  it  is  in  use.  The 
nearest  scattering  surface  (other  than  the  plywood  boxes)  is  the  desktop,  which  is 
18  inches  below  the  source-detector  plane.  Two  concrete  walls  and  a  metal  door  are  also 
quite  close  to  the  source,  as  shown  in  Figure  15.  This  high-scattering  geometry  is  not 
optimum  for  actually  performing  neutron  spectrometiy,  but  for  original  equipment  setup 
and  testing,  the  arrangement  is  adequate.  The  configuration  is  convenient  to  use  during 
construction  and  testing  because  it  can  be  easily  varied  on  the  desktop  as  improvements 
are  made  in  the  system.  For  radiation  safety,  the  operator's  station  is  located  in  another 
room,  approximately  4  meters  from  the  source. 
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Electronics:  Set-up  and  Optimization 

The  setup  and  optimization  of  the  electronics  in  the  spectrometry  system  is 
accomplished  in  two  parts.  In  the  first  part,  the  detection,  PSD,  and  output  subsystems 
are  connected  together  and  optimized  for  maximum  pulse  shape  discrimination.  This 
configuration  is  shown  in  Figure  16.  In  the  second  part  of  the  setup  process,  the 
electronics  are  reconnected  in  the  configuration  shown  previously  in  Figure  14.  The  two 
parts  of  the  setup  process  are  discussed  sequentially. 


Figure  16.  Electronics  configuration  for  PSD  work 


Pulse  Shape  Discrimination  (PSDV  With  the  electronics  configured  as  shown  in 
Figure  16,  the  HV,  preamp,  DLA,  and  PSA  are  each  optimized  in  turn.  Optimization 
consists  of  observing  the  PSD  picture  on  the  MCA,  shown  in  Figure  17,  and  adjusting  the 
controls  on  each  module  to  achieve  the  largest  separation  between  the  gamma  and  neutron 
peaks  and  to  achieve  the  best  resolution  of  the  peaks  as  measured  by  the  peak-to-valley 
ratios.  The  ADC  and  MCA  are  set  to  1024  channels  full  scale  for  best  resolution,  and  the 
ADC  gate  is  not  used.  The  linear  output  from  the  PSA  is  used  as  the  signal  input  to  the 
ADC. 


Figure  17.  Experimentally-determined  optimum  PSD  for  this  spectrometry  system. 
The  narrow  peak  at  the  left  corresponds  to  gamma  events,  while  the  broader  peak 
at  the  right  corresponds  to  neutron  events 
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The  HV  supply  must  be  set  to  negative  high  voltage  to  correspond  to  the  input 
voltage  requirements  of  the  PMT  and  PMT  base.  Keeping  in  mind  that  the  maximum  high 
voltage  allowed  for  the  tube  is  -2700  V,  the  HV  is  varied  in  increasing  steps  of  100  V, 
beginning  at  -1000  V.  It  is  found  that  the  PSD  picture  of  Figure  17  is  essentially 
impossible  to  see  at  -1000  V.  It  becomes  progressively  better  from  -1 100  to  -1300  V  then 
begins  to  show  anomalies  at  higher  voltages.  The  optimum  setting  is  -1300  V,  as  this 
setting  achieved  the  best  combinatio  n  of  peak  separation  (14  channels)  and  peak-to-va!ley 
ratio:'  (y-valley:  4.9,  n-valley:  9.8)  Tbt  value  of  the  high  voltage  is  read  from  the  digital 
display  on  the  HV  supply  rather  than  from  the  dial  settings,  under  the  assumption  that  the 
voltmeter  is  more  accurate  than  the  scale  on  the  dials.  The  discrepancy  between  the  two 
values  is  35  V,  with  the  display  reading  lower  than  the  dials. 

The  only  control  that  can  be  varied  on  the  preamplifier  is  the  input  capacitance, 
which  can  be  set  to  either  0,  100,  200,  500,  or  1000  pf  Observing  the  PSD  spectrum  on 
the  MCA  for  each  of  these  settings,  it  is  found  that  the  spectrum  becomes  progressively 
worse  with  increasing  capacitance,  until  the  PSD  spectrum  is  unintelligible  at  1000  pf. 

The  best  resolution,  then,  is  at  0  pf  (y-n  sep:  14  ch,  y-valley  ratio:  5.3  n-valley  ratio:  10). 
The  optimized  output  from  the  preamp,  as  viewed  on  the  oscilloscope,  looks  like  the  plots 
shown  in  Figures  18  and  19.  Figure  18  shows  a  collection  of  gamma  and  neutron  event 
pulses  of  different  amplitudes,  while  Figure  19  is  a  magnification  of  just  one  of  these 
positive,  unipolar  pulses  in  which  the  detail  of  the  pulse  can  be  seen.  These  plots  are 
direct  screen  dumps  from  the  oscilloscope. 
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Figure  1 8.  Oscilloscope  trace  of  output  from  the  preamp  showing  several  signals  from 
the  PuBe  UT022  y/n  source  (vertical  scale:  5  mV/div,  horizontal  scale:  1  ms/div) 

The  next  module  in  line  is  the  delay  line  amplifier.  Several  controls  must  be 
adjusted  on  the  DLA.  They  are  the  input  signal  polarity,  integration  time,  DC  offset,  pole 
zero,  and  gain.  The  input  polarity  switch  must  be  set  to  positive  to  match  the  positive 
input  signal  from  the  preamp.  An  integration  time  of  0.25  ^s  is  absolutely  necessary  to 
preserve  the  timing  information.  Setting  the  integration  time  to  either  0. 1  or  0.04  psec 
distorts  the  DLA  output  into  a  square  wave  and  results  in  no  discernible  PSD  on  the 
MCA.  The  correct  setting  of  the  DC  offset,  a  ten-turn  screw  adjustment,  is  verified  on  the 
oscilloscope  by  setting  the  0-V  line  of  the  DLA  output  to  be  the  same  as  the  grounded 


Figure  19,  Oscilloscope  trace  of  output  from  preamp  showing  single  signal  from 
PuBe  #T022  with  magnified  horizontal  scale  to  show  detail  of  pulse 
(vertical  scale:  5  mV/div,  horizontal  scale:  50  ps/div) 


oscilloscope  trace.  The  pole  zero,  another  ten-turn  screw  adjustment,  must  also  be 
adjusted  by  using  the  oscilloscope.  The  correct  setting  of  the  pole  zero  eliminates  over- 
and  under-shoot  of  the  0-V  line  at  the  end  of  the  DLA  output  pulse.  Finally,  the  DLA 
gain  is  adjusted  by  optimizing  the  PSD  spectrum  on  the  MCA.  While  PSD  ij  discernible 
over  a  range  of  gains  from  3  to  250,  the  optimum  setting  is  found  to  be  15  (y-n  sep:  14  ch, 
y-valley  ratio:  5.0,  n-valley  ratio;  9.9).  The  oscilloscope  trace  of  the  optimized  bipolar 
output  from  the  DLA  is  shown  in  Figure  20. 
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Figure  20.  Oscilloscope  trace  of  bipolar  output  from  DLA  showing  single  signal  from 
PuBe  #T022  (vertical  scale:  O.I  V/div,  horizontal  scale.  0.5  ;.  s/div) 


The  next  module,  which  is  the  PSA,  has  settings  for  time  range,  input  voltage 
discrimination,  walk  adjustment,  and  upper  and  lower  level  discriminators  for  the  time 
window.  The  first  thing  to  do  is  to  set  the  full-scale  time  range  from  to  0.4  ps.  This  both 
expands  the  PSD  picture  and  moves  it  toward  the  center  of  the  MCA  screen.  Increasing 
the  input  discrimination  makes  the  valley  between  the  gamma  and  neutron  peaks  more 
prominent.  Up  to  now,  the  input  discrimination  has  been  set  at  0.2  V,  but  by  raising  this 
control  to  0.7  V,  the  background  and  dead  time  are  both  lowered. 

The  walk  adjustment  of  the  PSA  is  very  important.  This  control  is  a  ten-turn 
screw  with  no  scale  on  the  module  to  indicate  its  position.  It  must  be  adjusted  by 
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watching  the  PSD  spectaim  on  the  MCA.  The  optimum  adjustment  of  this  control  is 
reached  when  the  gamma  and  neutron  peaks  reach  their  maximum  separation 
(y-n  sep:  17  ch,  y-valley  ratio:  70,  n-valley  ratio:  31).  Turning  the  screw  in  either 
direction  from  this  optimum  point  will  bring  the  peaks  closer  together,  or  even  smear  them 
together  into  a  single  peak.  If  no  PSD  is  observed  during  the  earlier  stages  of  this 
optimization  procedure,  try  adjusting  this  screw  and  then  repeating  the  procedure  with  the 
previous  modules.  Oscilloscope  traces  of  the  linear  output  of  the  PSA,  which  is  used  as 
the  signal  input  to  the  ADC  when  determining  PSD,  are  shown  in  Figures  21  and  22. 
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Figure  21  shows  a  single  pulse,  while  Figure  22  uses  a  five-times  compression  feature  of 
the  oscilloscope  to  show  five  consecutive  pulses  on  a  compressed  horizontal  scale.  The 
first  four  pulses  of  lower  amplitude  are  gamma  pulses,  while  the  last,  higher-amplitude 
pulse  corresponds  to  a  neutron  event. 
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Figure  22.  Oscilloscope  trace  of  linear  output  from  PSA  using  5x  compression  feature 
on  scope  to  show  4  gamma  pulses  (-4.3  V)  and  1  neutron  pulse  (-7.6  V) 
(vertical  scale:  I  V/div,  horizontal  scale;  2.5  ps/div)  1 


Finally,  the  neutron/gamma  discrimination  window  on  the  PSA  must  be  set. 


To  do 


this,  the  linear  output  from  the  PSA  is  used  as  the  signal  to  the  ADC,  while  the  in-window 
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logic  output  is  used  as  the  gate  input  to  the  ADC.  To  perform  this  setup,  the  PSA,  ADC, 
GDG,  and  delay  amplifier  must  all  be  adjusted  as  part  of  the  same  process.  It  is  necessary 
to  send  the  linear  output  of  t’'e  PSA  through  the  delay  amplifier  before  connecting  it  to  the 
ADC,  and  the  logic  signal  must  be  routed  through  the  GDG  before  it  is  sent  to  the  ADC. 

First,  the  GDG  must  be  set  so  that  the  logic  pulse  from  the  PSA  is  of  the  proper 
amplitude  and  duration  to  act  as  a  gate  for  the  ADC.  The  ADC  used  in  this  system 
requires  that  the  amplitude  of  the  gate  be  between  +3  and  +10  V  and  that  the  duration  be 
a  minimum  of  1  |as.  Further,  the  gate  input  must  precede  the  triggering  of  the  ADC's 
lower  level  discriminator  (LLD)  by  the  signal  pulse  by  approximately  100  ns  and  remain 
until  after  the  LLD  has  been  triggered  [25: 1-6].  The  GDG  has  three  controls,  the  delay, 
amplitude,  and  width  of  the  gate  pulse.  Set  the  delay  to  0. 1  ps,  which  is  the  minimum 
delay  available.  The  amplitude  and  width  controls  must  both  be  set  using  the  oscilloscope, 
since  they  have  no  marked  scales  on  the  face  of  the  module.  The  square  pulse  of  the  gate 
signal  is  set  to  3.4  V,  which  is  within  the  range  of  acceptable  amplitudes.  This  low  value 
is  chosen  because  there  is  a  large  ringing  at  the  leading  edge  of  the  gate  pulse.  By 
choosing  this  amplitude,  the  maximum  amplitude  of  the  leading  edge  ring  is  6.4  V,  which 
is  also  within  the  acceptable  range  of  amplitudes.  Therefore,  the  amplitude  is  in  range  no 
matter  which  part  of  the  gate  pulse  is  triggering  the  aDC.  The  width  of  the  gate  pulse  is 
set  at  3.0  V  in  order  to  comfortably  meet  the  timing  requirements  for  the  gate. 

Next,  the  delay  setting  on  the  delay  amplifier  must  be  set  so  that  the  PSD  signal 
pulse  falls  within  the  gate  signal.  To  set  this  control,  use  the  oscilloscope  and  add  delay  to 
the  signal  until  its  leading  edge  falls  at  least  100  ns  after  the  leading  edge  of  the  gate  pulse. 
A  setting  of  0.25  ps  places  the  leading  edge  of  the  signal  1  ps  after  the  leading  edge  of  the 
gate,  while  the  gate  continues  for  0.5  ps  after  the  signal  pulse  has  ended.  The  oscilloscope 
trace  showing  the  properly-set  PSD  signal  and  gate  pulses  is  presented  in  Figure  23. 
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Figure  23.  Oscilloscope  trace  showing  the  inputs  to  the  ADC  used  for  setting  PSD. 
Large  square  pulse  is  the  gate,  while  the  smaller  pulse  centered  in  the  gate 
is  the  PSD  signal  (vertical  scale;  1  V/div,  horizontal  scale:  1  ps/div) 


Now  that  the  gate  has  been  established,  the  discriminators  on  the  PSA  can  be  set. 
The  discriminators  will  be  set  so  that  gamma  events  fall  within  the  window,  while  neutron 
events  fall  above  the  window.  The  lower  level  discriminator  is  set  to  0  V  so  that  all 
low-amplitude  pulses  are  counted  as  gamma  events.  Now,  using  the  PSD  spectrum  on  the 
MCA,  set  the  coincidence  control  on  the  ADC  to  "coincidence"  and  raise  the  upper  level 
discriminator  of  the  PSA  from  zero  until  the  entire  gamma  peak  and  about  half  of  the 
valley  is  displayed  on  the  MCA.  To  verify  this  setting,  change  the  ADC  coincidence 
switch  to  "anti-coincidence"  and  look  at  the  PSD  spectrum  on  the  MCA.  Now  only  the 
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neutron  peak  and  the  other  half  of  the  valley  should  be  displayed.  The  optimum  setting  for 
the  ULD  on  the  PSA  is  4.34  V.  Looking  at  the  oscilloscope  trace  of  the  series  of  linear 
PSA  output  pulses  in  Figure  22,  one  can  see  that,  to  within  the  uncertainty  of  the  trace 
calibration  and  the  dial  setting,  this  discriminator  level  falls  above  the  gamma  pulses  and 
below  the  neutron  pulse. 

Energy.  The  second  part  of  the  electronics  setup,  which  involves  the  energy 
subsystem  and  part  of  the  output  subsystem,  is  much  less  complicated  than  the  work 
required  previously.  To  begin  this  process,  the  electronics  must  be  configured  back  to  the 
normal  spectrometry  setup  shown  in  Figure  14.  The  tasks  for  this  part  of  the  setup  are  to 
set  the  amplifier,  delay  amplifier,  ADC,  and  MCA  such  that  the  energy  signal  properly  falls 
within  the  timing  gate  and  so  that  the  energy  spectrum  meets  the  requirements  for  input  to 
the  FORIST  unfolding  code. 

First,  set  the  input  voltage  discriminator  on  the  PSA  to  its  lowest  setting,  which  is 
0.010  V,  in  order  to  record  the  low-amplitude  end  of  the  spectrum.  Next,  the  amplifier 
needs  to  be  adjusted  to  produce  a  good  Gaussian  pulse  for  the  energy  signal.  To  make 
these  adjustments,  view  the  output  of  the  amplifier  on  the  oscilloscope.  Now,  turn  off 
both  the  base  line  reject  and  delay  features,  since  they  are  not  needed  for  this  application. 
Then  set  the  shaping  time  of  the  pulse  to  0.5  fis.  Finally,  set  the  gain  to  approximately 
750.  This  gain  will  be  adjusted  later  during  the  energy  calibration  of  the  spectrometer. 

Next,  connect  the  output  of  the  GDG  to  the  other  input  of  the  oscilloscope  and 
view  the  PSD  gate  and  the  energy  signal  together  so  that  the  delay  of  the  energy  signal 
may  be  adjusted  on  the  delay  amp  to  meet  the  gating  requirements  of  the  ADC.  It  is 
found  that  setting  the  delay  amp  to  2.5  ps  meets  the  requirements.  The  oscilloscope  trace 
of  the  optimized  energy  signal  and  its  PSD  gate  are  shown  in  Figure  24. 
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Figure  24.  Oscilloscope  trace  of  inputs  to  the  ADC  for  recording  energy  spectra. 
The  large  square  pulse  is  the  gate  produced  by  the  PSD  subsystem,  while  the 
Gaussian  pulse  within  the  gate  is  the  signal  generated  by  the  energy  subsystem 
(vertical  scale;  1  V/div,  horizontal  scale;  1  ps/div) 


Finally,  the  number  of  channels  on  the  ADC  needs  to  be  set  to  256,  as  this  is  the 
number  of  channels  used  by  the  ZSHIFT  code  to  prepare  the  FORIST  input  file.  Now, 
raise  the  lower  level  discriminator  of  the  ADC  up  from  zero  until  the  dead  time  shown  on 
the  ADC  and  MCA  is  less  than  5%.  This  setting  is  found  to  be  at  0.22  V.  Set  the  upper 
level  discriminator  to  its  maximum  value,  which  is  10.04  V,  in  order  to  include  as  many  of 
the  pulses  as  possible  in  the  spectrum.  Now  set  the  zero  offset  control  to  approximately 
6.00.  The  position  of  this  control  will  be  finalized  during  the  energy  calibration  process. 
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Finally,  the  electronics  are  set  up  and  the  spectrometry  system  is  ready  for  use.  A 
summary  of  the  settings  for  all  of  the  modules  is  given  in  Table  3. 


Data  Collection  Procedure 

There  are  three  steps  to  be  performed  in  the  laboratory  each  time  a  data  run  is 
made  with  the  spectrometer.  These  are; 

1 .  verify  the  PSD  discrimination  point  set  on  the  PSA, 

2.  energy  calibrate  the  system,  and 

3.  take  the  four  spectral  data  sets  required  by  FORIST. 

Verifying  the  PSD  discrimination  point  is  not  a  difficult  task.  Simply  connect  the 
electronics  in  the  configuration  of  Figure  19  and  view  a  neutron/gamma  source  on  the 
MCA  with  the  coincidence  switch  on  the  ADC  set  to  "normal"  so  that  both  the  gamma 
and  the  neutron  peak  are  displayed.  The  flip  the  switch  to  "coincidence"  to  verify  that 
only  the  gamma  peak  is  displayed  and  flip  the  switch  to  "anti-coincidence"  to  see  that  only 
the  neutron  peak  is  displayed.  Throughout  the  two  months  that  tliis  spectrometer  was  in 
constant  operation,  no  drift  in  the  discrimination  point  was  seen  and  no  corrections  were 
needed. 

The  next  task  is  to  energy  calibrate  the  system  using  the  small  ^^Na  source  (#217). 
To  begin,  the  electronics  must  be  reconfigured  into  their  normal  setup  for  spectrometry,  as 
shown  in  Figure  14,  and  the  settings  on  the  modules  must  be  in  their  spectrometry 
positions,  as  listed  in  Table  3.  Set  the  coincidence  switch  on  the  ADC  to  "normal"  to 
disable  the  PSD  gate,  since  this  is  a  gamma-only  source.  Then  place  the  22Na  source  in 
contact  with  the  detector  to  achieve  the  maximum  counting  rate. 


52 


Table  3.  Summary  of  electronics  settings 


High  Voltage  Supply 

HV 

-1300  V 

Preamplifier 

Input  Capacitance 

0  pf 

Delay  Line  Amplifier 

Fine  Gain 

0.3 

Coarse  Gain 

50 

Integration  Time 

0.25  ns 

Input 

positive 

Pulse  Shape  Analyzer 

Time  Range 

0.4  (is 

Input  Discriminator 

0.700  V  for  PSD 

If  tl 

0.010  V  for  spectrometry 

Upper  Level  Discriminator 

4.34  V 

Lower  Level  Discriminator 

0.00  V 

Gate  and  Delay  Generator 

Delay 

0.10  |is 

Amplitude 

3.5  V 

Width 

3.0  |is 

Amplifier 

Fine  Gain 

0.750 

Coarse  Gain 

1000  (adjusted  during  runs) 

Shaping  Time 

0 . 5  }is 

Base  Line  Reject 

off 

Delay 

out 

Delay  Amplifier 

Delay 

1.0  ns  for  PSD 

It 

2.5  ns  for  spectrometry 

Analog-to-Digital  Converter 

Conversion  Gain 

1024  ch  for  PSD 

If  It 

256  ch  for  spectrometry 

Upper  Level  Discriminator 

10.04  V 

Lower  Level  Discriminator 

0.22  V 

Zero  Offset 

6.00 

Coincidence 

normal  for  full  PSD  spectrum 

tt 

coincidence  for  gammas 

It 

anticoincidence  for  neutrons 

The  energy  calibration  of  the  system  must  yield  four  results;  the  gain  calibrations 
of  the  spectrometer  at  both  the  high  and  low  gain  settings  required  for  FORIST  and  the 
zero  intercepts,  or  channel  locations  of  zero  energy,  for  both  gain  settings.  FORIST 
requires  the  gain  settings  to  be  near  or  slightly  above  0,00625  lu/ch  for  the  high  gain 
setting  and  0.0625  lu/ch  for  the  low  gain  setting,  while  also  requiring  both  of  the  zero 
intercept  locations  to  be  greater  than  zero. 

At  amplifier  gain  settings  above  approximately  150,  both  Compton  edges  from  the 
22Na  source  can  be  seen  on  the  MCA.  However,  at  amplifier  gain  settings  below  150,  the 
Compton  edge  of  the  lower-energy  gamma  (0.51 1  MeV)  is  no  longer  resolved  on  the 
MCA.  Figures  25  and  26  show  the  measured  Compton  recoil  electron  spectra  of  22Na 
#T028  at  the  high  (amp  gain  =  750)  and  low  (amp  gain  =  75)  gain  settings  used  by 
FORIST.  Because  the  Compton  edge  of  the  first  gamma  cannot  be  resolved  at  the  low 
gain  setting  required  by  FORIST,  a  linear  interpolation  and  least-squares  extrapolation 
computer  program,  called  CALIBER.  TK  has  been  written  to  calculate  all  of  the  gain  and 
zero  intercept  information  required  by  FORIST. 

To  use  CALIBER.  TK,  a  22Na  recoil  electron  spectrum  is  recorded  at  each  of  four 
different  amplifier  gain  settings,  two  of  which  correspond  to  the  high  and  low  gain 
requirements  of  FORIST  and  the  other  two  of  which  fall  in  between.  The  two 
intermediate  gain  runs  are  required  because  they,  along  with  the  high  gain  setting,  provide 
the  three  points  needed  for  the  least  ^squares  fit  which  will  be  extrapolated  to  the  low  gain 
setting  to  find  the  location  of  the  first  gamma's  Compton  edge  at  low  gain.  For  the 
spectrometry  system  constructed  for  |his  research,  the  following  four  amplifier  gain 
settings  are  used.  With  the  fine  gain  of  the  amplifier  to  0.750,  varying  the  coarse  gain  dial 
to  100,  200,  500,  and  1000  results  in  gains  of  75,  150,  325,  and  750  respectively.  The 
settings  at  75  and  750  are  a  factor  of  10  apart,  as  required  by  the  low  and  high  gain  inputs 
for  FORIST,  and  they  also  provide  the  necessary  gain  values  in  terms  of  lu/ch,  as  will  be 
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Figure  25.  Laboratory  measurement  of  ^^Na  #T028  Compton  electron  spectrum 


at  high  gain  setting 

shown  soon.  The  two  intermediate  settings  at  150  and  325  each  resolve  both  of  the  22Na 
Compton  edges,  so  they  can  be  used  in  the  least-squares  calculation. 

The  instructions  for  the  use  of  CALIBER.TK  are  divided  into  7  steps,  as  follows. 
The  complete  CALIBER. TK  model  is  included  in  Appendix  C. 

First,  record  the  four  spectra  described  above.  Use  ^^Na  source  #217,  placed 
directly  on  the  side  of  the  detector,  and  count  each  spectrum  for  2  minutes  to  achieve 
reasonable  counting  statistics.  Print  the  counts-per-channel  data  and  locate  on  these  print¬ 
outs  the  channel  locatio'  .s  of  the  peaks  corresponding  to  the  two  gammas 
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Figure  26.  Laboratory  measurement  of  ^^Na  #T028  Compton  electron  spectrum 


at  low  gain  setting 

(with  the  exception  of  the  first  gamma  peak  in  the  low  gain  data,  which  is  not  resolved). 

Second,  start  CALIBER. TK  and  display  the  list  sheet  on  the  screen.  Fill  th».  list 
called  "amp_gain"  with  the  four  settings  used  on  the  amplifier,  listing  them  from  lowest  to 
highest.  Then  fill  the  list  called  "peakl"  with  the  channel  locations  of  the  three  known 
locations  of  peak  1  corresponding  to  the  three  higher  gain  settings.  Fill  the  "peak2"  list 
with  the  channel  locations  of  the  second  peak  for  each  gain  setting.  Finally,  fill  the  two 
lists  called  "PkCtsl"  and  "PkCts2"  with  the  counts  recorded  in  each  of  the  known  peak 
locations. 
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Third,  display  the  rule  sheet  on  the  screen.  The  rules  are  divided  into  four 
sections,  which  will  be  used  one  at  a  time.  The  first  set  of  rules,  which  are  labeled  "find 
counts  in  1/2  height  of  peaks,"  compute  the  counts  in  the  half  height  location  of  each  of 
the  known  peaks  entered  in  the  previous  step.  To  use  the  first  set  of  rules,  make  sure  that 
there  is  an  asterisk  in  the  status  column  of  the  rule  sheet  next  to  the  two  rules  in  this 
section.  All  of  the  rules  in  the  other  sections  should  have  the  letter  "C"  next  to  them, 
indicating  that  they  are  currently  canceled  and  will  not  be  calculated.  Now  press  the  FIO 
key  to  solve  the  first  set  of  rules.  This  will  fill  the  two  lists  called  "Ctsl "  and  "Cts2"  with 
the  results  of  this  calculation: 

Fourth,  move  to  the  table  sheet  and  print  the  three  tables  titled  "peak_table", 
"pkl_cts_table",  and  "pk2_cts_table".  These  three  tables,  as  printed  by  CALIBER.TK  for 
one  of  the  calibration  runs  in  this  research,  are  shown  in  Figure  27.  Now  move  to  the  plot 
sheet  and  print  the  plot  titled  "peak_plot".  This  plot  serves  .s  a  visual  check  on  the  data 
used  in  the  calibration,  as  both  of  the  gamma  peak  position  vs.  amplifier  gain  lines  should 
be  linear.  If  the  plotted  points  do  not  appear  linear,  double-check  the  data  entered  at  the 
beginning  of  this  procedure  and  fix  any  errors  found.  The  plot  of  this  data  is  shown  in 
Figure  28. 

Fifth,  the  channel  locations  of  the  half-heights  of  the  Compton  edges  must  be 
computed  using  a  linear  interpolation  routine  and  the  results  of  the  half-height  calculation 
in  step  three.  To  begin  this  step,  move  to  the  rule  sheet.  Cancel  the  first  set  of  rules  and 
un-cancel  the  second  set  of  rules,  titled  "linear  interpolation  routine."  Now  set  the  screen 
so  the  variable  sheet  and  list  sheet  are  both  displayed.  In  the  variable  sheet,  find  the 
section  titled  "linear  interpolation  routine."  This  routine  must  be  run  seven  times,  once  for 
each  known  peak.  For  each  peak,  enter  the  counts  in  the  half-height,  calculated  in  step 
three,  on  the  input  line  of  the  variable  "y3".  Then  look  at  the  print-out  from  the  MCA  th.at 
contains  the  data  for  the  appropriate  peak  and  find  the  channels  with  counts  immediately 
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Peak  positions  vs  amp  gain 


1  Amp  gain 

peak  1  (cb 

peak  2  (ch)  | 

75 

17 

150 

12 

32 

375 

24 

70 

750 

43 

133 

Peak  1:  Counts  in  peak  and  1/2  height 


1  amp  gain 

1  peak  ch 

1  peak  cts 

1  1/2  Cts  1 

75 

150 

12 

27108 

13554 

375 

24 

12147 

6073.5 

750 

43 

6529 

3264.5 

Peak  2:  Counts  In  peak  and  1/2  height 

1  amp  gain 

peak  ch 

peak  cts 

1/2  Cts  1 

75 

17 

7602 

3801 

150 

32 

4012 

2006 

375 

70 

1711 

855.5 

750 

133 

923 

461.5 

Figure  27.  Three  tables  created  by  CALIBER, TK 


above  and  below  the  half-height  location.  Enter  the  channel  number  of  the  channel  to  the 
left  of  the  half-height  location  in  the  variable  "xl"  and  enter  the  counts  in  that  channel  in 
the  variable  "yl Then  enter  the  counts  in  the  channel  immediately  to  the  right  of  the 
half-height  location  in  the  variable  "y2".  Now  press  F9  and  the  linearly-interpolated 
location  of  the  half-height  will  appear  in  the  output  column  of  the  variable  "x3''.  Save  this 
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Figure  28.  Plot  of  peak  positions  vs.  amplifier  gain  to  verify  linearity 
of  the  data  in  Figure  27 


result  in  the  list  sheet  in  either  "pl_12"  or  "p2_12",  as  appropriate,  and  run  the  calculation 
again  until  the  half-height  locations  of  all  seven  peaks  have  been  computed  and  saved. 

Sixth,  display  both  the  variable  sheet  and  the  rule  sheet  on  the  screen.  In  the  rule 
sheet,  cancel  the  set  of  rules  and  un-cancel  the  third  set  of  rules,  titled  "computation  of 
peak  1  position  using  least  squares  fit  to  peak  1  data."  In  the  variable  sheet,  find  the 
section  titled  "peak  1  position  routine"  and  make  sur";  that  the  correct  number  of  spectra, 

4,  is  listed  at  the  variable  "n",  and  that  the  amplifier  setting  at  low  gain  (75  for  this  system) 
is  listed  at  the  variable  "x".  Now  press  F9  and  the  extrapolated  value  for  the  half-height  of 
the  first  peak  at  the  low  gain  setting  wili  appear  in  the  variable  sheet  in  the  output  column 


the  two  gammas  in  the  high  gain  setting  in  the  variables  "xl_h"  and  ’’x2_h",  and  enter  the 
half-height  positions  of  the  two  gammas  in  the  low  gain  setting  in  the  variables  "xl  1"  and 
"x2J".  Now  press  F9  and  the  slope,  y-intercept,  zero-intercept,  and  gain  (in  lu/ch)  will  be 
calculated  for  the  high  and  low  gain  settings.  The  part  of  the  CALIBER.  TK  variable  sheet 
which  lists  the  answers  is  shown  in  Figure  30. 


St  Input - 

Naae - 

•  341 

El 

N«V 

9»»al  aax  coxpton  alactron  anergy 

1  062 

E2 

MaV 

ganaal  aax  coap'.on  alactron  anargy 

48.388608 

xl_h 

ch 

location  of  gaual  at  high  gain 

143.33239 

*2_h 

ch 

location  of  gaaaal  at  high  gain 

•-h 

.00768298  HeV/ch 

slopa  at  high  gai:- 

b_h 

-.0307687  MaV 

y-intarcapt  at  high  gain 

tjti 

4.0047917  ch 

*aro  Intarcept  at  high  gain  (RESULT) 

C-h 

.00640348  lu/ch 

high  gain  (RESULT) 

9.7S23036 

Xl_l 

Ch 

location  of  gaaaal  at  low  gain 

18.864485 

x2_l 

ch 

location  of  gaaaal  at  low  gain 

■_l 

.07913398  MaV/ch 

slope  at  low  gain 

-.4306312  MaV 

y-intarcapt  at  low  gain 

5.4429117  ch 

xaro  intarcept  at  low  gain  (RESULT) 

C.l 

.06593665  lu/ch 

low  gain  (RESULT) 

Figure  30.  Part  of  CALIBER.TK  variable  sheet  showing  gain  and  zero-intercept  results 


FORIST  requires  that  the  resulting  zero-intercepts  must  both  be  greater  than  zero, 
and  that  the  gains  must  be  >  0.0625  lu/ch  and  Gff  >  0.00625  !u/ch.  These  conditions 
may  not  be  satisfied  with  the  first  calibration  attempt,  in  which  case  the  amplifier  gains 
must  be  adjusted  and  the  procedure  repeated.  Initially,  this  requires  some  trial-and-error, 
but  once  a  successful  calibration  is  made  all  future  calibrations  are  much  easier. 
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With  the  calibration  of  the  spectrometry  system  completed,  taking  data  is  quite 
simple.  FORIST  requires  that  a  foreground  (no  shadow  shield)  and  background  (with 
shadow  shield)  data  run  be  taken  at  the  calibrated  liigh  and  low  gain  settings,  for  a  total  of 
four  data  runs.  The  four  required  Compton  electron  spectra  of  the  PuBe  #M-1 1 70  source 
are  shown  in  Figures  3 1  to  34,  in  the  following  order:  high  gain  foreground,  high  gain 
background,  low  gain  foreground,  and  low  gain  background.  The  foreground  and 
background  runs  at  each  gain  have  the  same  general  shape,  but  there  are  many  more 
counts  in  the  foreground  spectra. 

The  counting  time  for  these  four  runs  must  be  long  enough  to  ensure  good 
counting  statistics.  For  the  sources  used  in  this  research,  these  times  were  30  minutes  to 
2  hours  for  the  small  PuBe  (#T022)  source,  15  minutes  for  the  large  PuBe  (#M-1 170) 
source,  and  1  hour  for  the  large  22Na  (#T028)  source.  The  four  spectral  data  sets  taken 
during  a  spectrometry  data  run  must  each  be  saved  to  disk  and  then  converted  to  ASCII 
format  using  the  "LOTCNV"  command  in  the  AccuSpec  MCA  software. 


CJ  O  3 


Figure  33.  PuBe  #M'l  170  Compton  recoil  electron  low  gain  foreground  spectrum 


2200 


Pu6c  tM1170,  run  8,  lou  gain  gammaa.  rou  dots,  IS  Nov  92 


Figure  34.  PuBe  #M-1 170  Compton  recoil  electron  low  gain  background  spectium 


Data  Processing  and  Unfolding  Procedure 

Now  that  the  data  has  been  collected  in  the  laboratory,  it  must  be  processed  and 
unfolded  to  yield  the  spectrum  of  the  incident  radiation.  The  steps  in  this  process  are; 

1 .  prepare  the  laboratory  data  for  input  to  ZSHIFT, 

2.  prepare  and  run  ZSHIFT, 

3.  run  FORIST,  and 

4.  plot  the  results  from  FORIST. 


The  first  step,  preparing  the  laboratory  data  for  input  to  ZSHIFT,  requires  first 
transferring  the  data  files  from  the  floppy  disk  where  they  are  stored  in  the  laboratory  to 
the  mainfi'ame  account  that  contains  ZSHIFT  and  FORIST.  The  next  requirement  is  to 
modify  the  four  ASCII-format,  recoil-spectra  data  files  so  that  they  may  be  read  using 
standard,  formatted  input  statements  in  FORTRAN.  The  data  files  contain  several  lines  of 
unneeded  information  before  the  spectral  data,  which  ZSHIFT  ignores.  The  spectral  data, 
then,  is  originally  formatted  into  rows  of  8  data  values,  with  each  data  value  separated 
from  its  neighbors  by  2  blank  spaces.  Since  the  data  values  contain  different  numbers  of 
digits  depending  on  their  magnitudes,  the  effect  of  this  format  is  to  give  the  data  a  random 
look  on  the  screen.  To  be  read  by  ZSHIFT,  these  four  spectra  must  be  modified  by  hand 
so  that  the  data  fall  into  8  columns  of  uniform  width  and  spacing,  a  process  which  takes 
about  twenty  minutes. 

The  second  step,  preparing  and  running  ZSHIFT,  involves  the  following  steps: 

1 .  identify  the  four  spectral  data  files, 

2.  enter  the  gain,  zero-intercept,  and  counting  time, 

3.  set  the  input  format  for  the  MCA  data  files, 

4.  choose  the  appropriate  response  matrix,  and 

5.  set  the  desired  parameters  for  FORIST. 

Referring  to  the  ZSHIFT  computer  code,  contained  in  Appendix  D,  all  required  user- 
inputs  and  brief  instructions  for  their  use  are  contained  on  the  first  two  pages  of  the  code. 

The  first  set  of  inputs,  which  identifies  the  spectral  data  files  to  be  unfolded, 
requires  that  the  user  give  four  sets  of  information.  First,  the  actual  filenames  of  the  four 
data  files  must  be  entered  as  the  values  "FILE(l)”  through  "FILE(4)".  These  entries  are 
used  to  find  and  open  the  four  data  files.  Second,  short,  6-character  titles  for  each  of  the 
four  files  must  be  entered  in  the  variables  "TITLE(l)"  through  "TITLE(4)".  These  titles 
are  used  to  label  the  four  spectra  in  some  of  the  output  generated  during  the  data 
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processing  and  unfolding  procedure.  Third,  longer  identifiers  of  up  to  80  characters  must 
be  entered  as  "!DENT(1)"  through  "IDENT(4)''.  These  descriptive  identifiers  are  also 
used  to  label  the  four  spectra  in  some  of  the  output  generated  during  the  data  processing 
and  unfolding  procedure.  Finally,  a  number  to  identify  the  laboratory  data  set  is  entered  in 
the  variable  "IDENTY".  This  integer  is  used  to  identify  the  data  in  some  of  the  output. 

The  second  set  of  inputs  gives  the  parameters  for  the  experimental  data.  The 
values  "ZI(1)"  and  "ZI(2)"  contain  the  zero-intercepts  of  the  high  and  low  gain  runs, 
respectively.  The  variables  "HGAIN"  and  "LGAIN"  contain  the  gains  of  the  high  gain  and 
low  gain  spectra,  in  units  of  light  units  per  channel.  Finally  the  four  variables 
"HFNORM",  "HENORM",  "LFNORM”,  and  "LBNORM"  contain  the  counting  time,  in 
seconds,  for  each  of  the  four  spectra.  Counting  times  for  each  set  of  four  spectral  data 
runs  in  this  reseaich  were  generally  of  equal  length,  but  they  can  be  of  different  lengths  if, 
for  instance,  counting  statistics  on  one  or  two  of  the  spectra  need  to  be  improved. 

The  third  group  of  inputs  sets  the  input  format  for  the  MCA  data  files.  The 
variable  "M"  is  the  number  of  channels  in  the  MCA.  The  variable  "SKIP"  is  a  format  line 
that  tells  ZSHIFT  how  many  lines  of  unnecessary  header  information  to  ignore  before 
beginning  to  read  the  spectral  data.  The  variable  "INFORM"  is  a  format  line  that  tells 
ZSHIFT  the  format  of  the  spectral  data  in  the  for  data  files.  For  instance, 

"INFORM  =  '(1X,8F7.0)"'  indicates  that  there  are  8  data  values  on  each  line  of  data,  each 
of  which  are  7  digits  (including  spaces)  wide.  It  is  important  to  verify  that  this  format  line 
matches  the  hand-modified  format  of  the  data  in  the  four  files. 

The  fourth  input  selects  the  response  matrix  that  FORIST  will  use  in  the  unfolding 
calculation.  Three  response  matrices  are  currently  available  for  FORIST.  The  "University 
of  Illinois  Gamma-Ray  Matrix"  is  contained  in  the  file  "GAMMA.MAT"  and  should  be 
used  for  all  gamma  unfolding  runs.  The  "Illinois  Neutron  Response  Matrix",  stored  in  the 
file  "NEUTRON.MAT"  was  not  used  during  this  research.  Instead,  the  "Illinois  Neutron 
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Response  Matrix  Modified  Below  2.5  MeV  by  NEH",  stored  in  the  file  "UINEUT.MAT" 
was  used  to  unfold  all  of  the  neutron  data.  This  choice  was  made  because  the  low-energy 
modifications  in  "UINEUT.MAT"  were  made  on  the  basis  of  more  recent  experimental 
data  and  computer  modeling. 

The  fifth  and  final  set  of  inputs  to  ZSHIFT  allows  the  user  to  select  the  desired 
parameters  for  the  FORJST  unfolding  routine.  These  selectable  parameters  are  "NBITS", 
"STRTCH",  "NPUNCH",  "ITRATE",  "ERRER",  and  "TAW".  The  "NBITS"  parameter 
is  the  number  of  bits  of  precision  of  the  computer  on  which  FORI  ST  is  running,  so  is 
therefore  machine-dependent.  The  values  for  several  machines  are  listed  in  the  ZSHIFT 
code.  I’he  "STRTCH"  parameter  is  a  multiplier  on  the  Gaussian  window  widths  used  in 
the  unfolding  routine.  This  parameter  should  be  near  1 .0  for  best  results.  The 
"NPUNCH"  parameter  tells  FORJST  which  unfolded  spectra  to  include  in  the  output  files. 
For  this  work,  "NPUNCH"  was  always  set  to  3,  so  that  the  smoothed  and  unsmoothed 
spectra  could  both  be  plotted  and  compared  against  one  another.  The  "ITRATE" 
parameter  sets  the  number  of  smoothing  iterations  performed  by  FORJST.  One 
smoothing  iteration  is  usually  sufficient  [16:340].  The  "ERRER"  parameter  allows  the 
user  to  choose  the  desired  statistical  error,  in  peicent,  of  the  unfolded  result.  The  value  of 
this  parameter  is  generally  set  at  3%.  [16:341]  The  final  parameter,  "TAW",  sets  the 
constraint  weight,  which  requires  the  unfolded  spectrum  to  be  "near  zero",  thus  helping  to 
dampen  out  oscillations.  This  parameter  works  best  when  set  between  1.0  and  10.0 
[14:35,82-83].  During  the  course  of  this  research,  the  parameters  "STRTCH", 

"ITRATE",  "ERRER",  and  "TAW"  were  systematically  varied  in  order  to  study  their 
effect  on  the  unfolded  spectra.  The  complete  set  of  FORIST  parameters  are  summarized 
in  Table  4. 
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Table  4.  FORIST  parameters  for  neutron  and  gamma  unfolding 


For  neutron  unfolding: 

NBITS 

23 

STRTCH 

1.0 

NPUNCH 

3 

ITRATE 

1 

ERRER 

3 

TAW 

10.0 

For  gamma  unfolding: 

NBITS 

23 

STRTCH 

1.2 

NPUNCH 

3 

ITRATE 

1 

ERRER 

3 

TAW 

1.0 

After  the  user  has  selected  all  of  the  ZSHIFT  input  and  run  ZSHIFT,  two  output 
files  are  generated.  The  first  of  these,  called  ZSHIFT.OUT,  summarizes  all  of  the  input 
and  output  spectral  data.  This  file  is  used  for  debugging  ZSHIFT  and  is  not  used  in  the 
actual  unfolding  process.  The  second  output  file,  called  SPECTRA.  IN,  is  the  input  file  for 
FORIST.  Samples  of  these  two  files,  corresponding  to  the  inputs  contained  in  the 
ZSHIFT  code  in  Appendix  D,  are  located  in  Appendices  G  and  H. 

After  the  spectral  data  has  been  prepared  by  ZSHIFT,  it  is  very  simple  to  unfold  it 
using  FORIST.  The  FORIST  code  itself  has  no  user  inputs,  but  simply  reads  the 
SPECTRA. IN  file  created  by  ZSHIFT,  along  with  the  appropriate  response  matrix 
specified  in  SPECTRA.IN.  FORIST  generates  three  output  files.  The  first  of  these,  called 
FORIST.OUT,  contains  a  great  deal  of  information  which  is  useful  in  debugging  the  code. 
The  second  output  file,  called  PUNCH.OUT,  is  intended  as  an  input  file  to  a  graphics 
routine  which  was  not  used  in  this  research,  so  this  file  is  not  used  in  the  current 
configuration  The  third  file,  called  TOGRAPH.OUT,  is  the  input  file  to  the  graphics 
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routine,  called  GRAPHICS. TK,  which  is  used  in  this  research.  Samples  of  FORIST.OUT 
and  TOGRAPH.OUT,  corresponding  to  the  ZSHIFT  files,  are  contained  in 
Appendices  I  and  J. 

The  final  step  in  the  data  processing  and  unfolding  process  is  to  graph  the  data  in 
the  TOGRAPH.OUT  file  using  the  GRAPHICS. TK  code.  To  do  this,  the 
TOGRAPH.OUT  file  must  first  be  transferred  from  the  mainframe  account  to  a  floppy 
disk,  since  the  GRAPHICS. TK  code  runs  on  a  personal  computer.  To  use 
GRAPHICS. TK,  start  the  program  and  then  load  the  TOGRAPH.OUT  data  using  the 
ASCII  loading  feature  of  TK.  Now  press  the  FI  0  key  to  convert  the  energy  scale  of  the 
data  from  eV  to  MeV.  The  unfolded  spectra  can  now  be  viewed  by  looking  at  the  two 
plots  in  the  model:  the  "unsmooth"  plot  contains  the  FORIST  results  before  the 
smoothing  routine  is  performed,  while  the  "smooth"  plot  contains  the  FORIST  results 
after  one  smoothing  iteration.  The  plots  should  be  titled  and  scaled  appropriately  before 
printing.  If  "NPUNCH"  or  "ITER"  are  changed  from  the  recommended  values  of  3  and  1, 
GRAPHICS. TK  will  have  to  be  modified  accordingly. 
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V.  Experimental  Results 


Now  that  the  theory  and  operation  of  the  spectrometry  system  have  been 
explained,  the  last  step  is  to  present  the  laboratory  results  used  to  validate  the  system. 
These  results  fall  into  four  major  groups:  the  test  of  the  modified  FORIST  code  with 
previously-validated  DT  neutron  data,  the  neutron  and  gamma  spectra  from  the  small 
PuBe  source  (#T022),  the  neutron  and  gamma  spectra  from  the  large  PuBe  source 
(#M1 170),  and  the  gamma  spectrum  from  the  large  22Na  source  (#T028). 


Test  of  Modified  FORIST  Code 

« 

Upon  modifying  the  FORIST  unfolding  code  for  use  on  the  VAX  6420  mainframe 
and  before  using  the  new  version  to  unfold  data  taken  with  the  newly-constructed 
spectrometry  system,  the  proper  functioning  of  the  FORIST  code  itself  needs  to  be 
verified.  To  make  this  verification,  a  recoil  proton  da*a  set  from  a  DT  neutron  generator 
was  provided  by  Hartley  [12;-].  This  set  consists  of  the  four  spectral  data  files  required  by 
FORIST.  Additionally,  the  appropriate  inputs  to  ZSHIFT,  summarized  in  Table  5,  were 
supplied  with  the  data.  The  unfolding  run  yielded  the  neutron  spectrum  shown  in  Figure 
35.  The  expected  DT  neutron  spectrum  for  this  data  set,  which  was  taken  at  a 
135°  angle  with  respect  to  the  beam,  is  a  single  peak  at  13.6  MeV  [12;-].  which  is  what  is 
obser/ed  in  the  plot.  The  x-axis  of  the  plot  shows  the  neutron  energy  in  MeV,  while  the 
y-axis  shows  the  relative  intensity  of  each  data  point,  in  units  of  flux  per  MeV.  A  solid 
line  connects  the  data  points,  while  vertical  error  bars  are  plotted  to  show  the  +/- 
uncertainty  in  the  relative  intensities. 
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Table  5,  Inputs  to  ZSHIFT  for  FORIST  verification  run 


Spectral  Parameters: 

HGAIN 

0.0064 

lu/ch 

LGAIN 

0.0720 

lu/ch 

ZI ( 1 )  (high  gain) 

2.37 

ch 

ZI (2)  (low  gain) 

3.05 

ch 

HFNORM 

298800 

HBNORM 

460000 

LFNORM 

192000 

LBNORM 

593000 

FORIST  Parameters: 

NBITS 

23 

STRTCH 

1.2 

NPUNCH 

3 

ITRATE 

1 

ERRER 

3 

TAW 

1.0 

Figure  35.  DT  neutron  spectrum  used  to  verify  correct  operation 
of  modified  FORIST  code 
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PuBe  Source  #T022 

The  recoil  proton  and  Compton  electron  spectra  of  PuBe  source  #T022  were 
measured  in  the  laboratory  using  the  specirometry  system  constructed  for  this  project. 
They  were  then  unfolded  to  yield  the  incident  neutron  and  gamma  spectra  from  the  source. 

The  ZSHIFT  input  parameters  for  the  neutron  data  set  are  given  in  Table  6  and  the 
resulting  unfolded  PuBe  neutron  spectrum  is  shown  in  Figure  36.  The  peaks  of  the 
experimentally-determined  spectrum  in  Figure  36  are  compared  with  the  accepted  values 
for  the  PuBe  peak  energies  [  1 :23 1  ]  in  Table  7.  The  agreement  between  the  measured 
spectrum  and  the  true  spectrum  is  within  the  5  -  10%  uncertainty  expected  with  the 
system.  The  unresolved  peak  at  1 .4  MeV  is  not  always  observed  in  the  accepted  PuBe 
spectrum  [1 :341],  while  the  observed  peak  at  5.3  MeV  is  probably  an  oscillation  in  the 
FORIST  unfolding  process.  The  plot  in  Figure  36  is  the  result  of  the  0**’  smoothing 
iteration  of  FORIST,  meaning  that  the  spectrum  has  been  smoothed  but  has  not  been 
iterated  on  using  FORIST's  "Iterated  Smoothing  Technique."  The  result  generated  after 
the  first  smoothing  iteration  was  not  used  because  FORIST  interpreted  most  of  the  peaks 
as  oscillations  and  tried  to  smooth  them  away. 


Table  6.  Input..,  ;o  ZSHIFT  for  PuBe  #T022  neutron  run 


HGAIN 

0.00631 

lu/ch 

LGAIN 

0.0650 

lu/ch 

ZI f 1 )  (high  gain) 

4.515 

ch 

ZI (2)  (low  gain) 

5.580 

ch 

HFNORM 

3600 

sec 

HBNORM 

3600 

sec 

LFNORM 

3600 

sec 

LBNORM 

3600 

sec 
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ruB#  ltT022  Neutron  Spectrum,  run  E,  18  Nov  92 


2  3 


N«utron  Ericrou  (M*V) 


Figure  36.  Experimentally-deten.iined  PuBe  #T022  neutron  spectrum 


Table  7.  Comparison  of  true  PuBe  neutron  enertries  to  measured  PuBe  #T022  spectrum 


True  Energy  (MeV)  Observed  Energy  (MeV)  Error  {%) 


1,4 

— 

— 

2.0 

2.0 

±0.0 

3.1 

3.1 

±0.0 

CO 

4.4 

-8.3 

— 

5.3 

— 

6.6 

’  6.2 

-6.1 

7.7 

7.2 

-6.5 

9.8 

9.6 

-2.0 

The  ZSHIFT  input  parameters  for  the  gamma  data  set  are  given  in  Table  8,  while 
the  resulting  unfolded  PuBe  gamma  spectrum  is  shown  in  Figure  37.  The  gamma  peak  of 
the  experimental  spectrum  is  at  4.25  MeV,  while  the  true  location  of  the  peak  is  at 
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4.43  MeV  [13:559],  for  an  error  of -4. 1%,  which  is  less  than  the  system  uncertainty  of 
5-10%.  The  results  of  the  first  smoothing  iteration  were  used,  as  they  gave  the  best 
resolution  of  the  peak. 


Table  8.  Inputs  to  ZSHIFT  for  PuBe  #T022  gamma  run 


HGAIN 

0.00631 

lu/ch 

LGAIN 

0.0650 

lu/ch 

ZI(1)  (high  gain) 

4.515 

ch 

ZI (2)  (low  gain) 

5.580 

ch 

HFNORM 

1800 

sec 

HBNORM 

1800 

sec 

LFNORM 

1800 

sec 

LBNORM 

1800 

sec 

Figure  37.  Experimentally-determined  PuBe  #T022  gam  ;a  spectrum 
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PuBe  Source  #M-1 170 

The  recoil  proton  and  Compton  electron  spectra  of  PuBe  source  #M-1 170  were 
also  measured  in  the  laboratory  and  unfolded  using  FORIST. 

The  ZSHIFT  input  parameters  for  the  neutron  data  set  are  given  in  Table  9  and  the 
resulting  unfolded  PuBe  neutron  spectrum  is  shown  in  Figure  38.  The  peaks  of  the 
experimentally-determined  spectrum  in  Figure  38  are  compared  with  the  accepted  values 
for  the  PuBe  peak  energies  [1 :23 1]  in  Table  10. 


Table  9.  Inputs  to  ZSHIFT  for  PuBe  #M-1 1 70  neutron  run 


HGAIN 

0.00640 

lu/ch 

LGAIN 

0.0663 

lu/ch 

ZI (1 )  (high  gain) 

4.005 

ch 

ZI (2)  (low  gain) 

5.509 

ch 

HFNORM 

900 

sec 

HBNORM 

900 

sec 

LFNORM 

900 

sec 

LBNORM 

900 

sec 

The  agreement  between  the  true  and  measured  spectra  is,  again,  less  than  the  / 

5-10%  uncertainty  associated  with  the  spectrometry  system.  The  peak  at  1 .4  MeV 
remains  unresolved,  while  peaks  caused  by  FORIST  oscillations  show  up  at  4.2  and  7.2 
MeV,  The  plot  in  Figure  38  is  the  0^  smoothing  iteration  plot,  as  the  peak  information 
was  smoothed  out  in  the  first  smoothing  iteration.  Note  from  the  error  bars  that  even 
though  the  counting  time  for  this  spectrum  was  just  one  fourth  of  the  counting  time  for  the 
#T022  neutron  spectrum,  the  error  was  smaller  for  the  #M-1 170  source.  This  is  because 
the  #M- 1 1 70  source  was  enough  stronger  that  several  times  more  counts  were  recorded  in 
the  lab,  even  though  the  counthig  time  was  much  shorter. 
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Figure  38.  Experimentally-determined  PuBe  #M-n70  neutron  spectrum 


Table  10.  Comparison  of  true  PuBe  neutron  energies  to  measured 

_ PuBe  #M-1 170  spectrum _ 

I  True  Energy  (MeV)  Observed  Energy  (MevT  Error  (TT 


1.4 

2.0  2.0  +0.0 

3.1  3.0  -3.2 

4.2 

4.8  4.9  +2.1 

6.6  6.4  -3.0 

7.2 

7.7  7.8  +1.3 

9.8  10.0  +2.0 


The  ZSHIFT  input  parameters  for  the  gamma  data  set  are  given  in  Table  1 1,  while 
the  resulting  unfolded  gamma  spectrum  is  shown  in  Figure  39.  The  gamma  peak  of  the 


78 


experimental  spectrum  is  at  4.40  MeV,  while  the  true  location  of  the  peak  is  at  4.43  ivleV 
[13:559],  for  an  error  of  -0.7%,  which  is  again  well  within  the  uncertainty  limits  of  the 
spectrometry  system.  The  features  below  3.5  MeV  are  FORJST  oscillations.  Again 
notice  that  the  error  bounds  on  the  spectrum  are  much  smaller  than  for  the  #T022  result, 
due  to  the  better  counting  statistics.  The  results  of  the  first  smoothing  were  again  used  for 
the  gamma  plot. 


Figure  39.  Experimentally-determined  PuBe  #M-1 170  gamma  spectrum 
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Table  1 1 .  Inputs  to  ZSHIFT  for  PuBe  #M-1 1 70  gamma  run 


HGAIN 

0.00631 

lu/ch 

LGAIN 

C.0650 

lu/ch 

ZI ( 1 )  (high  gain) 

4.515 

ch 

ZI (2)  (low  gain) 

5.580 

ch 

HFNORM 

1800 

sec 

HBNORM 

1800 

sec 

LFNORM 

1800 

sec 

LBNORM 

1800 

sec 

Na-22  Source  #T028 


As  a  final  check  on  the  operation  of  the  spectrometry  system,  the  gamma  spectrum 
from  the  large  22Na  source  was  measured  and  unfolded.  The  ZSHIFT  input  parameters 
for  the  run  are  given  in  Table  12,  while  the  resulting  unfolded  22Na  gamma  spectrum  is 
shown  in  Figure  40.  The  0.51 1  MeV  peak  was  not  resolved  because  the  FORIST  gamma 
response  matrix  does  not  go  below  0.8  MeV,  but  the  higher-energy  peak  was  very-well 
resolved.  The  measured  energy  for  the  second  peak  was  1.275  MeV,  which  was  exactly 
the  energy  expected.  The  counting  statistics  were  also  very  good  for  this  run,  as  is  shown 
by  the  very  small  error  bars.  The  plot  in  Figure  40  is  from  the  first  smoothing  iteration  of 
FORIST,  as  were  the  other  gamma  results  presented  previously. 


Table  12.  Inputs  to  ZSHIFT  for  ^^Na  #T02S  gamma  run 


HGAIN 

1 

0.00646 

lu/ch 

LGAIN 

0.0655 

lu/ch 

ZI  (1) 

(high 

gain 

i)  0.926 

ch 

ZI  (2) 

( low 

gain) 

2.180 

ch 

HFNORM 

3600 

sec 

HBNOPM 

i  3600 

sec 

LFNORM 

3600 

sec 

LBNORM 

3600 

sec 
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Figure  40.  Experimentally-determuied  #T028  gamma  spectrum 
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VI.  Summary  and  Conclusions 


In  summary,  the  goal  of  this  thesis  was  to  construct  and  test  a  neutron  and  gamma 
spectrometry  system  and  to  write  a  detailed  manual  for  the  use  of  this  system.  The 
detector  used  in  the  system  was  a  xylene-based  organic  scintillator  and  the  distinction 
between  gamma  and  neutron  events  in  the  scintillator  was  made  on  the  basis  of  pulse 
shape  discrimination.  The  recoil  proton  and  Compton  recoil  electron  spectra  recorded 
with  this  system  were  unfolded  using  the  FORIST  computer  code. 

This  thesis  first  discussed  the  theory  of  pulse  shape  discrimination,  spectrum 
unfolding,  and  energy  calibration.  It  then  described  all  of  the  apparatus  used  in  the 
spectrometry  system,  including  all  of  the  computer  software  associated  with  the  system. 
Then  a  detailed  user’s  guide  to  the  system  was  presented,  from  initial  setup  through 
plotting  the  unfolded  results.  Finally,  the  results  of  the  tests  made  on  the  system  were 
given. 

The  experimental  results  can  be  summarized  in  the  following  way.  The  unfolding 
of  Hartley's  DT  data  set  verified  that  FORIST  still  produced  the  correct  results  after  the 
modifications  made  to  the  code  in  order  to  allow  it  to  run  on  the  VAX  6420  mainframe 
The  results  of  the  three  gamma  data  runs  (PuBe  #T022,  PuBe  #M-1 170,  and  #T028) 
indicate  that  the  spectrometry  system  is  correctly  detecting  and  unfolding  incident  gamma 
events.  The  results  of  the  two  neutron  data  runs  (PuBe  #T022  and  PuBe  #M-1 170)  also 
show  good  agreement  with  the  accepted  PuBe  neutron  spectrum,  so  the  system  is 
correctly  detecting  and  unfolding  incident  neutron  events,  as  well.  The  discrepancies  in 
shape  between  the  measured  spectra  and  the  accepted  spectrum  are,  in  all  likelihood,  due 
to  the  high-scattering  environment  in  which  the  measurements  were  made. 


The  high  scattering  geometry  would  be  much  more  likely  to  affect  the  neutron  spectrum 
than  the  gamma  spectrum. 

Some  further  tests  with  other  neutron  sources  are  recommended  as  new  sources 
become  available.  The  PuBe  neutron  spectrum  is  quite  complicated,  so  it  is  not  the  easiest 
source  to  use  as  an  initial  test  for  a  neutron  spectrometer.  Better  sources  to  use,  because 
of  their  relative  simplicity,  would  be  a  DT  neutron  generator  and  a  252cf  spontaneous- 
fission  source. 

A  better,  lower-scattering  geometry  also  should  be  used  in  future  runs  with  this 
spectrometer.  The  room  in  which  the  thesis  measurements  were  made  was  good  from  a 
standpoint  of  radiation  safety,  but  the  proximity  of  the  desktop  and  walls  made  it  far  from 
optimum  from  a  low-scattering  point  of  view.  A  better  setup  would  be  to  put  the  source, 
detector,  and  shadow  shield  on  tripods  in  the  center  of  a  much  larger  room,  in  order  to  get 
as  far  as  possible  from  all  scattering  sources. 

Future  work  should  also  be  done  with  the  FORIST  and  ZSHIFT  codes.  As  it 
currently  exists  at  the  end  of  this  thesis,  the  codes  run  successfully  on  the  VAX  6420 
mainframe  and  produce  an  output  file  which  can  be  plotted  using  GRAPHICS. TK.  The 
codes  have  also  been  commented  to  make  them  ea'^ier  to  understand.  However,  it  is  the 
opinion  of  this  author  that  much  more  work  could  be  done  with  these  codes,  such  as 
improving  the  comments,  improving  the  logic  structure,  and  improving  the  plots  and  data 
generated  in  the  FORIST.OUT  file.  Especially  tedious  in  the  current  setup  is  the  need  to 
reformat  the  four  input  spectra  by  hand  before  running  ZSHIFT.  Ultimately,  it  would  be 
much  more  convenient  if  the  FORTRAN  programs  could  be  made  to  run  on  a  PC.  This 
would  allow  all  of  the  FORTRAN  and  TK  data  processing  software  to  be  located  on  the 
same  machine  that  is  used  to  record  the  MCA  data  in  the  lab.  This  work  could  even  turn 
out  to  be  the  subject  of  a  future  thesis. 


Finally,  now  that  this  spectrometry  system  has  been  constructed  and  tested,  the  Air 
Force  Nuclear  Engineering  Center  at  the  Air  Force  Institute  of  Technology  has  a  new  tool 
to  use  in  future  research  projects.  These  projects  could  include  such  things  as  building 
and  calibrating  an  associated  particle  system  for  the  Center's  DT  neutron  generator  or 
designing  the  Center's  new  hot  cell  facility  to  handle  a  strong  252cf  source.  This 
spectrometry  system  should  turn  out  to  be  a  tool  that  will  serve  the  Center  in  many 
different  applications  for  years  to  come. 
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Appendix  A:  Equipment  List 


The  following  is  a  list  of  the  equipment  used  to  build  the  detector  and  electronic 
circuit  for  the  spectrometry  system 


Table  13  Equipment  list 


ITEM 

Liquid  Scintillator 

5.0  cm  diameter 

5.0  cm  thickneos 

y  '  r  ^ 

ID  NUMBER 
5-9599-04 

Plexiglass  Light  Fire 
5.0  cm  diam.eter 

1.6  cm.  thickness 

Am  J r. -- 1 

r  *  r: 

none 

Photomultiplier  Tube 

Ha.m.am.atsu 

R32?-C2 

HC4”65 

PMT  base 

Crtec 

265 

1959 

High  Voltage  Supply 

Crtec 

556 

1928 

Preamplifier 

Ortec 

113 

6508 

Delay  Line  Amplifier 

Ortec 

460 

2014 

Pulse  Shape  Analyzer 

Ortec 

458 

567 

Gate  and  Delay 
Generator 

Ortec 

416A 

3639 

Delay  Amplifier 

Ortec 

427A 

1272 

Amplifier 

Ortec 

572 

47]  6 

Analog- to-Digital 
Converter 

Nuclear  Data 

ND570 

76-1234 

Multi-Channel 

Analyzer 

Nuclear  Data 

AccuSpec 

Ver  01 

Oscilloscope 

LeCroy 

9410 

2435 
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Appendix  B:  The  Plutonium-Berviliurr  Source 


Two  2'’9piutonium-^Beryllium  (a,n)  sources  were  used  in  tiiis  research  project  to 
provide  neutrons  in  the  presence  of  gammas.  The  first  was  a  small,  78-mCi  source 
(Monsar.to  Research  Corporation,  #T022)  in  a  tantalum  and  stainless  steel  cylinder  with  a 
diameter  of  1 .25  cm  and  a  height  of  1 .40  cm.  The  cylinder  was  attached  to  the  end  of  a 
plastic  rod  for  handling,  and  the  detector  was  placed  along  the  axis  of  symmetry  of  the 
cylinder  opposite  to  the  rod.  The  second  was  a  larger,  4.7-Ci  source  (Monsanto  Research 
Corporation,  #M-1 170)  in  a  tantalum  and  stainless  steel  cylinder  with  a  diameter  of 
2.59  cm  and  a  height  of  14.02  cm  [22: 1  ].  The  cylinder  was  attached  to  the  end  of  an 
aluminum  rod  for  handling,  and  the  detector  was  placed  perpendicular  to  the  axis  of 
symmetry  of  the  cylinder. 

The  stable  isotope  of  beryllium,  ^Be,  has  a  relatively  loosely-bound  neutron 
(1.7  MeV  binding  energy).  If  a  typical  alpha  particle  from  a  radioactive  decay  (5-6  MeV) 
strikes  a  ^Be  nucleus,  a  neutron  can  be  released  : 

+  %e  ->  +  n 

If  the  beryllium  is  mixed  with  a  long-lived  alpha-emitting  material,  such  as  239pu,  there 
will  be  a  nearly  constant  rate  of  neutron  production  [19:445]. 

However,  the  neutron  spectrum  from  this  source  is  not  a  simple  one. 
Plutonium-239  can  alpha  decay  by  emitting  an  alpha  particle  with  any  one  of  several 
different  energies.  The  primary  energies  and  probabilities  are  5. 1554  MeV  (73.3%), 
5.1429  MeV  (15.1%),  and  5.1046  MeV  (1 1.5%),  but  the  range  extends  from  4.640  to 
5.1554  MeV  [21:1464].  Additionally,  the  product  to  which  239pu  alpha-decays  (U-235) 
also  decays  by  alpha  emission.  A  process  of  alpha  and  beta  decay  continues  until  a  stable 
isotope  of  lead  is  reached.  To  further  complicate  the  alpha  emission  spectrum,  small 


amounts  of  contaminant  alpha  activity,  present  in  either  the  original  sample,  produced 
through  the  decay  of  a  precursor,  or  produced  through  decay  of  a  nuclide  produced  in  a 
neutron  capture  reaction  can  influence  the  overall  neutron  yield.  The  isotope  2'**Pu  is 
particularly  significant,  because  it  beta  decays  with  a  half-life  of  13.2  years  to  form 
Amencium-241,  which  is  a  strong  alpha-emitter.  The  neutron  yield  of  this  source  can 
therefore  gradually  increase  with  time  as  the  Am  accumulates  in  the  source.  An  original 
Pu-241  isotopic  fraction  of  0.7%  will  result  in  a  growth  rate  of  the  neutron  yield  of  2% 
per  year  until  an  equilibrium  state  is  reached  [17:24].  The  radioactive  production/decay 
series  is  shown  in  Figure  4 1 . 
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Figure  41 .  Plutonium  decay  series  [2:367] 
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From  2’’9py  i(5  products,  there  are  alphas  emitted  within  an  approxim.ate 

energy  range  of  1  to  6  MeV.  The  ^Be(a,n)*2c  reaaion  can  leave  the  product  carbon 
nucleus  in  its  ground  state  or  in  its  first  or  second  excited  states.  These  three  states  of  the 
product  carbon  atom  result  in  reaction  (^-values  of  5.704,  1 .29,  and  -1 .95  MeV, 
respectively  [30:436].  Therefore,  we  expect  to  find  neutrons  with  an  energy  spectrum 
from  0  to  almost  12  MeV  (£q  +  Q).  The  neutrons  produced,  then,  are  not  monoenergetic 

for  several  reasons.  These  reasons  are  [19:445]: 

1 .  the  many  diflferent  energies  of  the  alphas  produced, 

2.  the  slowing  of  alphas  that  will  occur  by  collision  in  the  source  material, 

3.  the  various  directions  of  neutron  emission  that  can  occur  relative  to  the 
incoming  alphas  (whose  directions  are  not  known), 

4.  the  possibility  that  the  is  left  in  one  of  three  different  energy  states,  and 

5.  the  possibility  that  the  ^Be  participates  in  a  multibody  breakup  reaction; 

'•a  +  %e  ->*a  +  *Be  +  n 

which  produces  neutrons  of  a  different  energy  than  those  from  a  two-body 
breakup  [30:435]. 

Furthermore,  once  the  neutrons  are  produced,  their  energies  may  be  modified 
before  leaving  the  source  by  one  or  more  of  the  following  processes  [20: 175]; 

1.  elastic  scattering  with  materials  in  the  source, 

2.  the  ^Be(n,2n)  reaction: 

n  +  ^Be  *Be  +  n  +  n 

3.  or  the  23^Pu(r^/)  fission  reaction: 

n  +  239pu  fission  fragments  +  (1  to  4)  n 
Vijaya  and  Kumar  [30:-]  and  Kumar  and  Nagarajan  [20;-]  calculated  the 
theoretical  neutron  spectrum  from  a  ^Be(a,n)  source  in  the  following  manner.  The 
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neutron  energy  distributions  from  the  ^Bc(a,n)>2c  reaction,  leaving  the  in  different 
energy  states,  were  calculated  using  the  equation  [30:436]; 


(t'D, 


1  1  da  ^ 


dE„ 


(14) 


where 

i  =  index  of  energy  states  of 

Ngg  =  number  density  of  beryllium  atoms  in  the  source, 

=  lowest  value  of  alpha  energy, 

=  highest  value  of  alpha  energy, 

£(£'^  =  stopping  power  for  alpha  particles  of  energy  E^^  in  the  source  material, 
{</<7/i/f2(0,£„)}_  =  the  center-of-mass  differential  cross-section  at  alpha  energy 

Ea  and  center-of-mass  angle  0  for  the  /***  excitation. 

The  term  /S£„(E^  is  given  by  [30:436]: 


iK^Mc) 


iK  +  Mc) 


(15) 


where  M„,  and  A/^  are  the  masses  of  the  neutron,  alpha  particle,  and  carbon  atom, 

respectively. 

The  results  of  this  (lalculation,/i(£„),^(£„),  and  f^{E„),  are  the  neutron  energy 
distributions  from  the  ^Be(lt,n)*2c  reaction  leaving  the  in  the  ground,  first,  and 
second-excited  states,  respectively.  Similarly,  the  neutron  energy  distribution  from  the 
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multibody  breakup  reaction,  which  gives  rise  to  a  low  energy  neutron  continuum,  is 
computed  using  [30:437], 


/« 


1 

da 

dE„ 

-{^:J 


dE^ 


(16) 


where  \doldE„{Efj^\£„  is  the  energy  differential  cross  section  at  bombarding  energy  E^ 
resulting  in  a  neutron  of  energy  E„.  The  primary  calculated  neutron  spectrum  is  then 
obtained  by  adding/1,/2,/3,  ^-nufcontinuum^  as  shown  in  Figure  42. 

Modifications  to  the  neutron  spectrum  caused  by  the  three  processes  listed  above 
were  made  to  the  primary  spectrum  with  the  help  of  a  Monte  Carlo  code.  Elastic 
scattering  causes  smearing  of  peaks,  while  the  (n,2n)  and  (n/)  reactions  boost  the  lov/ 
energy  part  of  the  spectrum. 

A  PuBe  neutron  source  has  a  neutron  production  rate  of  about  2x10^  neutrons  per 
second  for  each  Ci  of  239pu  [19:446],  so  the  /8-mCi  source  is  producing  approximately 
1 .6x105  neutrons  per  second,  while  the  4.7-Ci  source  is  producing  approximately 
9.4x10^  neutrons  per  second. 

In  addition  to  producing  neutrons,  the  PuBe  source  produces  a  significant  gamma 
ray  spectrum.  These  gammas  result  from  neutron  capture  reactions  in  the  heavy  nuclides 
in  the  source,  de-excitation  of  product  nuclides  after  an  alpha  decay,  and  energetic  scatters 
of  neutrons  with  materials  in  the  source  and  source  casing.  The  primary  gamma  energy  is 
at  4.43  MeV  and  is  caused  by  the  de-excitation  of  fi-om  its  first  excited  state  to  its 
ground  state  after  it  has  been  produced  in  the  ^Be((X,n)*2c  reaction  [13:559;21 :1 1], 
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■“  "  VARIABLE  SHEET  » 

St  Input -  Kaa« -  Output -  Unit - 


0 

PkCtal 

0 

PkCtal 

Ctal 

Ctal 

48 

xl 

ch 

3418 

yi 

cts 

X2 

ch 

3023 

ys 

Ctft 

X3 

ch 

3264.5 

y3 

Ct8 

a 

cts/ch 

4 

n 

75 

X 

A 

B 

y 

ch 

.341 

El 

M«V 

1.062 

E2 

M«V 

48.388608 

xl_h 

ch 

142.23239 

X2_h 

ch 

B_h 

.00768298 

M«v/ch 

b_h 

-.0307687 

M8V 

t_h 

4.0047917 

ch 

c_h 

.00640248 

lu/ch 

9.7522036 

xl_l 

ch 

18.864485 

x2_l 

ch 

■-1 

.07912398 

M«v/ch 

b_l 

-.4306332 

lf«V 

*_1 

5.4425117 

ch 

e_i 

.06593665 

lu/ch 

'  ' "  "  for  Acadeaic  Usa  Only 

-  Coaaent - 

Calibar.TX  7  Nov  92 

by  Robert  S.  Pope 
A  PROGRAN  TO  CALIBRATE  THE  PSD 
NEUTRON/CAIMA  SPECTROMETER  USING 
INFORMATION  FROM  THE  NA-22  SOURCE 

- CALCULATE  1/2  HEIGHT  OF  PEAK - 

eta  In  peak  1 
cte  in  peak  2 

cts  In  1/2  height  ol  peak  1 
eta  in  1/2  height  of  peak  2 

- LINEAR  INTERPOLATIOH  ROUTINE - 

high  point  channel 

high  point  counta 

low  point  channal  (calc  cs  xl  a  1) 

low  point  counts 

deairad  point  channal  (ANSWER) 

desired  point  eciuits 

alopa 

-—"PEAK  1  POSITION  ROUTINE - - 

nuBber  of  MCA  spectra 
input:  aap.gain  •  75 
output  froa  LaastSquara  function 
output  froa  LeastSguara  function 
daaired  result  —  peak  1  1/2  HT 

- GAIN  ROUTINE - 

gaaaal  aax  coapton  electron  energy 
gaaaa2  aax  coapton  electron  energy 

location  of  gaaaal  at  high  gain 
location  of  gaaaal  at  high  gain 
alope  at  high  gain 
y-intercept  at  high  gain 
sero  intercept  at  high  gain  (RESULT) 
high  gain  (RESULT) 

location  of  gaaaal  at  low  gain 
location  of  gaaaal  at  low  gain 
alopa  at  low  gain 
y-intercept  at  low  gain 
*ero  intercept  at  low  gain  (RESULT) 
low  gain  (RESULT) 


Naae - 

■■P-g*in 

peakl 

peak! 

PkCtal 

PkCtal 

Ctal 

Cta2 


pi-13 

p2_l2 


■leaenta —  Unit- 

4 
4 
4 
4 
4 

4 
4 
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•  Coaaent  ■  ■-  - . -  ■ 

—INPUT  THESE  S  LISTS  FROM  THE  22NA  DATA— 

ORTEC  572  Aaplifier  gain  netting 

channel  location  of  peak  1 

channel  location  of  peak  2 

counta  in  peak  1 

counts  in  peak  2 

—THESE  2  LISTS  ARE  OUTPUT  BY  1/2  HT  ROUTINE— 
counta  in  1/2  height  of  peak  1 
counta  in  1/2  height  of  peak  2 
—FILL  THESE  2  LISTS  BY  HAND  WITH 

RESULTS  OF  LINEAR  INTERPOLATION  ROUTINE— 
Channel  location  of  1/2  height  of  peek  1 

channel  location  of  1/2  height  of  peak  2 
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RULE  SHEET  ■  mi  . -  For  Acadealc  U«*  Only 

S  Rul*— . . . . — _ _ _ _ 1 

■--riHO  COUNTS  IN  1/J  HEIGHT  OP  PEAXS 

•  Cttl  -  PkCtNl  •  1/2 

•  CtaJ  ■  Pkct»2  •  1/2 

■—LINEAR  INTERPOLATION  ROUTINE 
C  ■  •  (X2  -  xl)  ■  y2  -  yl 
C  ■  •  (xl  -  xl)  -  y3  -  yl 
C  x2  >  xl  ♦  1 

■— COKPUTATION  OF  PEAR  1  POSITION  USING  LEAST  SQUARES  FIT  TO  PEAK  1  DATA 
C  Call  LeastSquare(n;A,B) 

C  y  -  A  ♦  B*x 

•—CAIN  ROUTINE 

C  N.h  -  (E2  -  El)  /  (X2_h  -  Xl_h) 

C  b_h  -  E2  -  ■.h  *  X2_h 
C  i_h  "  -b_h  /  ■_h 
C  CJ.  -  0.885  /  (x2_h  -  *_h) 

C  ■_!  -  (E2  -  El)  /  (X2_l  -  Xl_l) 

C  b_l  -  E2  -  ■_!  *  X2_l 
C  1_1  -  -b_l  /  ■_! 

C  0_1  «  0.885  /  (X2_l  •  1_1) 


——————  PROCEDURE  FUNCTION:  LaattSquara  m—  For  Acadaaic  Uaa  Only 

Couant:  Laast  Sqviara  function  to  calculata  paakl  1/2  loc.  at  C_1 

Paraaatar  Varlablaa: 

Input  Varlablaa:  n 

Output  Varlablaa;  A,B 

S  Stataaant - 

8U«X  ■  0 
auaY  ■  0 

au»X2  -  0 
aunXY  ■  0 
for  1  «  2  to  n 

avukX  -  auBX  *  'aap.qalnti] 

BURY  -  auaY  4  'pi.i2(l} 
auaX3  •  auaJ(2  *  (‘aap.qalrCl])'! 
auaXY  ■  auaXY  *  'aap.qalntl]  •  'pl_12Il] 
naxt  1 

count  ■  n  -  1 

dalta  •  count  *  au8X2  -  auRX‘2 
A  ■  (tuaX2  *  auaY  -  auaX  •  auaXY)  /  dalta 
B  •  (count  •  auaXY  •  auaX  •  auaY)  /  delta 


Naaa - 

paal(_plot 

qaaaaljplot 


PU3T  SHEET 


Plot  Type —  Olaplay  Option-  Output  Davlca 
Lina  chart  l.ECA 
Lina  chart  l.ECA 


For  Acadaalc  Uaa  Only 

Title - — 

Peak  poaitlona  va  aap 
Peak  1:  peak  and  1/2 


"  TABLE  SHEET  — 

Naaa -  Title - - - - 

pa*k_tabla  Peak  poaitlona  va  aap  qaln 

pkl_cta_tabla  Peak  It  Counta  In  peak  and  1/2  halqht 
pk2_ctB_tabla  Peak  2:  Counta  In  peak  and  1/2  halqht 
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Appendix  D:  ZSHIFT.FQR  Computer  Code 
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PROGRAM  Z SHI FT 


This  program  creates  an  input  file  for  the  FORIST  spectrum 
unfolding  code.  This  file,  called  SPECTRA. IN,  contains: 

-  the  name  of  the  response  matrix  file 

-  input  parameters  which  govern  the  operation  of  FORIST 

-  the  gain  and  zero  shift  for  each  spectrum 

-  the  four  compressed  and  shifted  spectra. 

Input  files  used  by  ZSHIFT.FOR; 

-  four  files  containing  the  MCA  spectral  data  from  the 
laboratory 

Output  files  created  by  ZSHIFT.FOR: 

-  ZSHIFT.OUT,  which  is  intended  for  debugging  purposes 

-  SPECTRA. IN,  which  is  the  input  file  for  FORIST 

User-inputs  required  by  ZSHIFT.FOR: 

-  all  user-inputs  are  in  the  two  boxed  sections  below. 

There  is  a  section  for  spectral  information  and 
identification  and  a  section  for  FORIST  parameters. 

Original  ZSHIFT.FOR  code  provided  by  Dr.  Nolan  E.  Hertel. 
Updated  and  modified  by  Robert  S.  Pope  in  November  1992. 


INTEGER  ERRER 
REAL  LFNORM, LBNORM, LGAiN 
CHARACTER* 6  TITLE 
CHARACTER* 12  FILE 

CHARACTER* 80  SKIP, INFORM, FMTl, FMT2, IDENT, RMATRIX 

DIMENSION  Y(513) ,YC(513) , IDENT ( 4 ) , TITLE (4 ) , FILE ( 4 ) , ZI ( 4 ) , SHC (4 ) 

OPEN  (UNIT=6,FILE='ZSHIFT.OUT'  , STATUS=’NEW* ) 

OPEN  (UNIT=7 , FILE= ' SPECTRA. IN ' , STATUS= 'NEW ' ) 


Spectrum  identification  and  information 
Input/Output  identification: 

FILE(n)  are  the  4  files  containing  the  MCA  spectral  data.  These 
must  be  in  the  order: 

1  =  high  gain  foreground 

2  =  low  gain  foreground 

3  =  high  gain  background 

4  =  low  gain  background 

TITLE (n)  are  the  6-character  labels  used  in  SPECTRA. IN  and 
FORI ST. OUT 

IDENT (n)  are  the  80-character  titles  used  in  SPECTRA. IN  and 
FORIST. OUT 

IDENTY  =  the  ID  number  of  the  data  run 

FILE(l)  =  ’PUBE_R9.HFN’ 

FILE (2)  =  ’ PUBE_R9 . LFN ' 

FILE (3)  =  ' PUBE_R9 . HBN ' 

FILE (4)  =  ' PUBE_R9 . LBN ' 

TITLE (1)  =  'PB9HFN' 
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TITLE(2)  =  'PB9LFN' 

TITLEO)  =  ’PB9HBN' 

TITLE (4)  =  •PB9LBN' 

IDENT(l)  =  'PuBe  #M1170  run  9  neutrons  hgf,  24  Nov  92' 

IDENT(2)  =  'PuBe  #M1170  run  9  neutrons  Igf,  24  Nov  92' 

IDENTO)  =  'PuBe  #M1170  run  9  neutrons  hgb,  24  Nov  92' 

IDENT(4)  =  'PuBe  #M1170  run  9  neutrons  Igb,  24  Nov  92' 

I DENT Y  =  9 

C  **  Zero  intercept  in  channels  [l=high  gain,  2=low  gain] : 

ZI (1)  =  4.005 
ZI(2)  =  5.509 

C  **  Gain  in  light  units  per  channel: 

HGAIN  =  0.00640 
LGAIN  =  0.0663 

C  **  Normalization  —  counting  time  in  seconds  for  each  spectrum: 
HFNORM  =  000900 
HBNORM  =  000900 
LFNORM  =  000900 
LBNORM  =  000900 

C  **  Format  of  data  to  be  input  from  MCA  data  files: 

C  SKIP  =  format  line,  one  "/"  for  each  header  line  to  be  skipped 

C  INFORM  =  format  line,  each  line  of  MCA  data 

C  M  =  Number  of  MCA  channels  in  input  spectra: 

C  ** 

SKIP  =  •  (/////////////////////////////////////////////////////)  ' 
INFORM  =  ' (1X,8F7.0) ' 

M  =  256 

C  **  Name  of  response  matrix  file  to  be  used  by  FORIST: 

C  UINEUT.MAT  or  NEUTRON. MAT  for  neutrons,  GAMMA. MAT  for  gammas 

C  ** 

RMATRIX  =  'UINEUT.MAT' 

C 

C  Parameters  for  FORIST 
C 

C  **  Number  of  bits  of  machine  precision: 

C  NBITS  =  36  for  CDC  1604 

C  NBITS  =  27  for  IBM  7090 

C  NBITS  =  23  for  VAX  6420 

C  NBITS  =  22  for  IBM  360 

C  ** 

NBITS=23 

C  **  Window  width  multiplier: 

STRTCH=1.0 

C  **  Punching  option:  governs  spectral  data  output  in  PUNCH. OUT 
C  and  TOGRAPH.OUT  files  from  FORIST 

C  0  =  none 

C  1  =  smoothed  spectrum  data 

C  2  =  unsmoothed  spectrum  data 

C  3  =  both  smoothed  and  unsmoothed 

C  *♦ 

NPUNCH=3 

C  **  Number  of  smoothing  iterations; 

ITRATE=1 

C  **  Desired  statistical  error  in  percent: 

ERRER=3 

C  **  Constraint  weight: 

TAW=10.0 

C - - - 
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FMTl  ='DYNMAT  1000  10  ( IX, A6, 14 , 10F7 . 0) ' 

FMT2  =’DYNMAT  2000  10  ( IX, A6, 14 , 10F7 . 0 ) ' 

ZI(3)  =  ZI(l) 

ZI (4)  =  ZI  (2) 

SHC(l)  =71(1)  -  0.5 
SHC(2)  =  ZI  (2)  -  0.5 
SHC(3)  =  SHC(l) 

SHC(4)  =  SHC(2) 

WRITE (6, 698)  RMATRIX 

698  F0R.MAT(1X,  'NAME  OF  RESPONSE  MATRIX  USED  IN  PROBLEM  =',A80) 
WRITE (7, 699)  RMATRIX 

699  FORMAT! IX, A80) 

WRITE (7, 700)  NBITS,STRTCH, IDENTY, NPUNCH, ITRATE, ERRER,TAW 

700  FORMAT(1X,I10,F10.3,4I10,F10.4' 

WRITE (7, 701)  HFNORM, LFNORM 
WRITE (7, 701)  HBNORM, LBNORM 
WRITE (7, 709)  HGAIN,LGAIN 
WRITE(7,709)  HGAIN,LGAIN 

701  FORM;^T(1X,2F10.0) 

709  FORMAT(1X,2F10.4) 

DO  100,  ITER=1,4 

OPEN  (UNIT=ITER, FILE=FILE ( ITER) , STATUS= ’ OLD ' ) 
READ(ITER,SKIP) 

READ (ITER, INFORM) (Y(I) , 1=1, M) 

WRITE (6, 601 )M, ZI (ITER) 

601  FORMAT (lOX, 13, IX, 'CHANNEL  ORIGINAL  SPECTRUM'/ 

+  10X,'ZERO  INTERCEPT  =',F7.3) 

WRITE ( 6,  620 ) TITLE ( ITER) 

620  FORMAT (1 OX, A6) 

WRITE(6,602) (Y(I),I=1,M) 

CALL  SHIFTY(Y,YC,M,SHC(ITER) ) 

WRITE(6,603)M 
WRITE(6,602) (YC (I) , I=1,M) 

602  FORMAT (8F10.0) 

603  FORMAT (lOX, 13, IX,  'CHANNEL  ZERO-SHIFTED  SPECTRUM') 

N=200 

IF (MOD (ITER, 2)  .NE.  0)  THEN 
WRITE(7,702)  FMTl 
ELSE 

WRITE(7,702)  FMT2 
END  IF 

WRITE (7, 702)  IDENT(ITER) 

702  FORMAT! IX, A80) 

DO  10  1=1, N 
Y(I)=YC(I) 

10  CONTINUE 
WRITE(6,604) 

604  FORMAT (lOX,  '200  CHANNEL  SUMMED  AND  SHIFTED  SPECTRUM' 

+  ,'  FOR  INPUT  TO  FORIST.FOR') 

NM=N-9 

DO  11  IB=1,NM, 10 
IE=IB+9 

WRITE(6,605)  TITLE ( ITER) , IB, (Y (I ) , I=IB, IE) 

605  FORMAT(1X,A6,1X,I3,1X,10F8.0) 

WRITE(7,606)TITLE(ITER) ,IB, ( Y ( I ) , I=IB, IE) 

606  FORMAT(lX,Ae,lX,I3,10F7.0) 

11  CONTINUE 
CLOSE  (ITER)  - 
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IF (MOD (ITER, 2)  .EQ.  0)  WR1TE(7,703) 
100  CONTINUE 

WRITE(7,703) 

703  FORMAT (IX) 

CLOSE  (6) 

CLOSE  (7) 

STOP 

END 


C  SUBROUTINE  SHIFTY 

C  This  routine  performs  a  zero  shift  on  a  pulse  height  spectrum 

C  using  a  linear  interpolation  scheme  borrowed  from  the  ALPHA-M 

C  code.  The  resulting  spectrum  is  then  compressed  into  200  channels 

C - 

SUBROUTINE  SHI FTY ( Y, YC, M, SHC) 

DIMENSION  X(513) ,Y(513) ,XC(513) ,YC{513) 

NP  =  0 

NP  =  INT(ABS(SHC) )  +  1 

C - if  zero  intercept  is  positive,  skip  the  following 


IF(SHC  .LT.  0.0) 
DO  10  1=1, M 

THEN 

YC(I+NP)  = 

Y(I) 

10 

CONTINUE 

DO  11  1=1, NP 
YC(I)=0.0 

11 

CONTINUE 

DO  12  1=1, M 
Y(I)=YC(I) 

12 

CONTINUE 

SHC=SHC+1.*N? 


END  IF 


DO  7  1=1, M 
XC(I)=I 
YC(I)=0.0 
X(I)=I*1.0  -  SHC 
7  CONTINUE 
DO  60  1=1, M 
DO  40  J=1,M 

IF(XC(I)-X(J) .EQ.O)  THEN 
45  YC(I)=Y(J) 

GO  TO  60 

ELSE  IF(XC(I)-X(J) .LT.O)  THEN 
50  YC(I)  =  (Y(J)  -Y(J-l)  )/{X(.-)-X(J-l)  ) 

YC(I)=Y(J-1)  +YC(I) * (XC(I)-X(J-l) ) 
GO  TO  60 
END  IF 

40  CONTINUE 
60  CONTINUE 
RFTURN 
END 


Appendix  E:  FORIST  FOR  Computer  Code 
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PROGRAM  FORIST 


FORIST  -  MAIN 

A  modification  of  the  COOLC  (FERDOR)  unfolding  code 
by  R.H.  Johnson,  Nuclear  Engineering  Program,  University  of 
Illinois,  13  Feb  1975. 

Modified  by  N.E.  Kertel  at  UT  Austin  to  read  data  related  to 
the  response  matrix  and  smoothing  windows  from  tape  UNIT  3.  Also 
modified  to  do  only  one  spectrum  per  run. 

Modified  once  more  by  N.E.  Hertel  on  2  Apr  1980  to  allow  TAW 
to  be  read  in  from  the  input  data  file  (UNIT  5) . 

Modified  by  R.S.  Pope  at  the  Air  Force  Institute  of 
Technology  in  October  1992  to  make  the  code  run  on  the  VAX  and  to 
make  the  code  FORTRAN  77  standard.  Also  reorganized  and  better- 
commented. 

PURPOSE: 

The  FORIST  code  is  used  to  unfold  data  from  an  NE-213  liquid 
scintillator  detector  in  order  to  generate  the  incident  gamma  or 
neutron  spectrum. 

INPUTS: 

All  user-selectable  inputs  are  contained  in  a  file  called 
SPECTRA. IN,  which  is  generated  by  the  ZSHIFT.FOR  program.  This 
file  contains  the  four  spectra  (high  gain  foreground,  low  gain 
foreground,  high  gain  background,  and  low  gain  background) 
measured  in  the  laboratory  and  the  identification  and  information 
about  these  spectra.  In  addition,  this  file  contains  all  user- 
selectable  parameters  for  FORIST,  as  well  as  the  filename  of  the 
appropriate  input  matrix. 

For  example,  if  FORIST  is  to  unfold  gammas,  the  two  input 
files  would  be: 

-  SPECTRA. IN  —  the  user  input  file 

-  GAi<MA.MAT  —  the  response  matrix  for  gammas 

Details  of  the  inputs  to  FORIST  can  be  found  in  the  ZSHIFT.FOR 
code. 

OUTPUTS : 

The  FORIST  code  generates  four  output  files.  The  first, 
WORKING. 10  is  not  intended  to  be  read  by  the  user.  It  is  simply 
a  file  used  during  matrix  manipulations  in  order  to  use  less 
memory  during  execution. 

The  second,  FORIST. OUT,  is  a  25  to  30  page  file  which 
contains  information  about  the  response  matrix  used  and  about  the 
input  spectra.  The  file  then  lists  many  of  the  results  calculated 
by  the  FORIST  program  during  exeo'.tion.  This  file  also  contains 
several  graphs,  none  of  which  seem  to  be  particularly  useful.  The 
best  use  for  this  file  is  in  debugging. 

The  third  file,  PUNCH. OUT,  contains  data  which  can  be  used  to 
plot  the  unfolded  spectrum  computed  by  FORIST.  This  particular 
file  is  formatted  to  be  read  by  a  plotting  routine  on  the  CYBER. 
Contact  Dr.  N.  E.  Hertel  for  information  about  this  plotting 
routine. 
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C  The  fourth  file,  TOGRAPH.OUT,  also  contains  data  which  can  be 

C  used  to  plot  the  unfolded  spectrum  computed  by  FOUIST.  This  file 

C  was  created  during  the  1992  modification  and  is  meant  to  be  read 
C  by  the  conimercially-available  program  "TK  SOLVER!".  This  file 

C  contains  the  energy  in  eV,  the  relative  counts 

C  for  each  energy  bin,  and  the  error  in  the  relative  counts.  This 

C  information  is  contained  in  single  columns.  The  spectral 

C  information  output  to  this  file  is  dependent  upon  the  punching 

C  option  chosen  in  the  program  ZSHIFT.FOR.  The  last  character  of 

C  each  variable  name  in  the  file  is  an  integer.  Like  integers  are 

C  grouped  together  and,  if  there  is  more  than  one  spectrum  in  the 

C  file,  the  highest-valued  number  is  the  least-smoothed  spectrum, 

C  continuing  in  decending  order  to  zero. 

C 

C  FOR  ADDITIONAL  INFORMATION: 

C  "A  User's  Manual  tor  COOLC  and  FORIST"  by  R.  H.  Johnson, 

C  Purdue  University,  IN  (1975)  is  available  through  RSIC  and  is 
C  considered  the  standard  manual  on  the  FORIST  code.  It  develops 

C  the  theory  nicely,  but  leaves  much  to  be  desired  in  the  area  of 

C  actually  operating  the  code. 

C  "Construction  and  Testing  of  A  Neutron  and  Gamma  Spectrometry 

C  System  Using  Pulse  Shape  Discrimination  with  an^  Organic 

C  Scintillator"  (AFIT/GNE/ENP/93M-6)  by  R.  S.  Pope,  Air  Force 

C  Institute  of  Technology,  Dayton,  OH  (1992)  is  a  thesis  which  more 
C  clearly  presents  the  actual  operation  of  the  FORIST  program  and 

C  associated  software  (ZSHFIT.FOR  and  TK  SOLVER). 

C 

C - - - 

DIMENSION  JDENT (20), IDENT (20) 

DIMENSION  A(156, 95), HT(156, 95) 

DIMENSION  B(156) ,S(156) ,UT(156) ,LLO(156) ,LUP(156) 

DIMENSION  VLO(156) ,VUP(156) ,FORTE(156) 

DIMENSION  W(95) ,Q(95) ,ATEMP(95) ,ELAB(95) ,X(95) 

DIMENSION  PLO(lll),PUP(lll),WINW(lll),PLAB(lll),ERSAVE(lll) 
DIMENSION  WIND (111), ERRF( 111) 

DIME.:SION  PLABLE(111,2) 

DIMENSION  BS{156,2) 

COMMON  A,HT 

INTEGER  ERRER  i 

CHARACTER* 80  RMATRIX 
C 

EQUIVALENCE (ATEMP, X) 

EQUIVALENCE ( PLABLE (1,1), PLAB ( 1 ) ) 

EQUIVALENCE ( PLABLE (1,2), WINW ( 1 ) ) 

EQUIVALENCE(BS(1,1),B(1) ) 

EQUIVALENCE ( BS ( 1 , 2 ) , S ( 1 ) ) 

C 

OPEN  (UNIT=2, FI LE=' WORKING. 10’  , STATUS= ’NEW ', FORM= ' UNFORMATTED ' ) 
OPEN  (UNIT=5,FILE=' SPECTRA. IN’  , STATUS= ’ OLD ’ ) 

OPEN  (UNIT=6,FILE=’FORIST.OUT’  , STATUS^ ’ NEW ’ ) 

OPEN  {UNIT=7,FILE=’ PUNCH. OUT’  , STATUS=’NEW’ ) 

OPEN  (UNIT=8, FILE=’TOGRAPH.OUT’,STATUS=’NEW’ ) 

READ (5, 113)  RMATRIX, 

113  FORMAT (IX, ABO) 

OPEN  (UNIT=3, FILE=  RMATRIX  , STATU3= ’ OLD ’ ) 

C 

KR=111 


I 
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MR  =  156 

READ  (3,905)  JDENT 
WRITE  (6,906)  JDENT 

READ (5,  903 ) MBITS,  STRTCH,  IDENTY, NPUNCH, ITRATE, ERRER, TAW 
READ(3, 901)NR,NC,NW, IL 

WRITE (6, 902 ) STRTCH, MBITS, NR, NC, NW, IDENTY, NPUNCH, ITRATE, ERRER, IL 

ITSA\/E=ITRATE 

IDENTY  =  IDENTY  -  1 

ERROR=ERRER*0.01 

CALL  MATIN ( FLAB, KR,NW, 2) 

C 

C  The  window  width  is  modified  by  the  factor  STRTCH. 

C 

DO  3  I  =  1,NW 

WINW(I)  =  STRTCH  *  WINW(I) 

WIND (I)  =  WINW(I) 

3  CONTINUE 

CAI,L  MATIN  (ELAB,NC,NC,1) 

CALL  AREAD(A,MR,NR,NC,NW) 

CALL  BSREAD (B, S, MR, NR, NC,  IDENT, LLO, LUP, VLO, VUP, FORTE, IL) 

ITRATE-- ITSAVE 
DO  5  JJ=1,NW 

WINW(JJ)  =  WIND(JJ) 

5  CONTINUE 

TDENTY  =  IDENTY  +  1 
C 

C  Calculate  Q  vector 

C 

DO  20  J=1,NC 
Q(J)=1.0E35 
DO  20  1=1, NR 

IF(A(I, J) )10,20,10 

10  QMIN=(AMAX1(B(I),0.0)+S(I))/A(I,J) 

IF(Q(J)-QMIN)20,20,15 
15  Q(J)=(2MIN 

20  CONTINUE 
C 

C  Calculate  HT  and  euclidean  norm  of  HT 
C 

ENORM  =0.0 
DO  25  1=1, NR 
DO  25  J=1,NC 

HT(I, J)=(Q(J)*A(I, J) )/S(I) 

ENORM  =  ENORM+HT(.I,  J)  **2 
25  CONTINUE 
WRITE  (2)  A 
BACKSPACE  2 

CALL  GINV(HT,MR,NR,NC,TAW/SQRT(FLOAT(NC) ) ,A,ATEMP) 

READ  (2)  A 
BACKSPACE  2 
DO  50  1=1, NR 
DO  50  J=1,NC 

HT(I,J)  =  (Q(J)*HT(I,  J) )/S(I) 

50  CONTINUE 

Calculate  the  dual  vector  D.  In  the  FERDOR  code,  it  is  not 
necessary  to  ma)<e  any  assumptions  about  the  so-called  true 
solution.  However,  there  exist  a  set  of  coefficients. 
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C  D{J),  J=1,NC  which  give  teh  variation  of  the  PLO  and  PUP  results 

C  as  the  W^k)  elements  of  the  window  matrix  are  varied.  These 

C  coefficients  are  sort  of  a  pseudo-solution. 

C  Then  calculate  the  weighted  sum  of  the  negative  dual  vector 

C  components.  In  ordinary  least  squares,  the  sum  of  squares  of 

C  residuals  is  minimized,  but  in  FERDOR  and  additional  term  1= 

C  implied  in  the  minimization  which  tends  to  prevent  the  pseudo- 

C  vector  from  becoming  negative.  If  it  does  go  negative  very  much, 

C  it  indicates  that  the  data  supplied  to  the  problem  are 

C  inconsistent  in  the  sense  that  no  non-negative  solution  exists. 

C  The  sum  of  negative  duals  indicates  the  degree  of  inconsistency. 

C  The  sum  should  be  between  0  and  NC  in  a  valid  problem. 

C 

ESL  =0. 

DO  110  J  =  l.NC 
X(J)  =  0. 

DO  104  I  =  1,NR 

X(J)  =  X(J)  +  HT(I,J)  *  B(I) 

104  CONTINUE 

SL  =  X(J)/Q(J) 

PCTSL  =  100.  *  SL 
ESL  =  ESL  +  SL**2 
110  CONTINUE 
C 

C  Calculate  BADJ,  the  adjusted  B  vector  which  is  consistent 

C  with  the  FERDOR  final  solution. 

C  Then  calculate  the  sum  of  squares  of  deviations  {B(I)- 

C  BADJ(I)),  1=1, NR  (weighted  by  1/S(I)).  This  sum  is  an  indicator 

C  of  the  validity  o  fthe  input  error  estimates,  S{I).  In  an 
C  ordinary  least  squares  problem,  the  usm  of  sqauares  of  deviations 

C  has  a  chi  square  distribution  which  should  approximately  be  equal 

C  to  NR.  In  the  FERDOR  formulation  (which  is  a  form  of  constrained 
C  least  squares),  it  is  ok  if  the  weighted  sum  of  squares  of 

C  deviations  is  less  than  NR,  but  it  indicates  inconsistency  if  it 

C  is  too  large.  The  factor  SQRT (weighted  sum  of  squares) /NR  is  the 

C  average  factor  by  which  the  input  errors  need  to  be  increased  in 

C  order  to  be  consistent. 

C 

WRITE  (6,906)  IDENT 
WRITE (6, 918)  TAW 
WRITE  (6,911) 

ESD  =  0. 

DO  120  I  =  1,NR 
BADJ  =  0. 

DO  114  J  =  1,NC 

BADJ  =  BADJ  +  A(I,J)  *  X(J) 

114  CONTINUE 

DEV  =  (B(I)  -  BADJ)/S(I) 

PCTDEV  =  100.  *  DEV 
ESD  =  ESD  +  DEV**2 

WRITE (6, 912) I, VLO( I), B (I)  ,S (I),  BADJ,  PCTDEV,  ESD 
120  CONTINUE 
FNC  =  NC 
FNR  =  NR 

DICAl  =  SQRT(ENORM) *FNC*0.5**NBITS 
DICA2  =  SQRT(ESL/FNC) 

DICA3  =  SQRT(ESD/FNR) 

SUMIND  =  DICA1+DICA2+DICA3 
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WRITE  (6,913)  DICA1,DICA2,DICA3,SUMIND 
WRITE (6, 918)  TAW 


Begin  window  loop. 

IF  (SUMIND  .GE.  5.0)  THEN 

58  WRITE  (6,909) 

END  IF 

59  WRITE  (6,906)  IDENT 
WRITE(6,918)  TAW 
WRITE  (6,907) 

IF  (NPUNCH.EQ. 1.0R.NPUNCH.EQ.3)  WRITE(8,920)  ' UTB ' , ITRATE, ' 

DO  500  K=1,NW 

CALL  WREAD(W,NC,NW,K,ELAB,PLAB,WINW) 

DO  60  1=1, NR 
UT(I)=0.0 


60 

CONTINUE 

DO  100  J=1,NC 

IF(W(  J))70,100,70 

70 

DO  80  1=1, NR 

UT(I)=UT(I)+W(  J)*HT(I,J) 

80 

CONTINUE 

100 

CONTINUE 

USSUM=0.0 

UTB=0.0 

DO  150  1=1, NR 

UTB=UTB+UT(I) *B(I) 

USSUM=USSUM+ (UT ( I ) *S ( I ) ) **2 
150  CONTINUE 

WUAQ=0.0 
DO  200  J=1,NC 
SUM=0 . 0 
DO  175  1=1, NR 

SUM=SUM+UT(I)*A(I, J) 

175  CONTINUE 

WUAQ=WUAQ+Q ( J) *ABS (W ( J) -SUM) 

200  CONTINUE 

EPSLON=SQRT (USSUM) +WUAQ 
PHIUP=UTB+EPSLON 
PHlLO=UTB-EPSLON 
ERR1=SQRT (USSUM) 

ERR2=WUAQ 

WRITE  (6,904)  PLAB (K) , K, PHILO, PHIUP, UTB, WI’JW (K) , ERRl, ERR2 
ERSAVE(K)  =  EPSLON 

I F (NPUNCH , EQ . 1 . OR . NPUNCH . EQ . 3 )  THEN 

WRITE ( 7, 914 ) PLAB ( K) , K, UTB, EPSLON, IDENTY, WINW ( K) 

IF  (K  .NE.  NW)  THEN 
WRITE (8, 921)  UTB,',’ 

ELSE 

WRITE (8, 922)  UTB 
END  IF 
END  IF 

PUP(K)  =  PHIUP 
PLO(K)  =  PHILO 
ERRF(K)  =  EPSLON  /  UTB 
500  CONTINUE 
C 

I F (NPUNCH . EQ . 1 . OR . NPUNCH . EQ . 3 )  THEN 
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\ 


WRITE (8, 9 '3) 

WRITE(8,920)  'EMTRATE,':' 
DO  600  I=1,NW-1 

WRITE(8,921)  PLAB(I),',' 
600  CONTINUE 

WRITE(8,922)  PLAB(NW) 
WRITE(8, 923) 

END  IF 
C 


IF(NPUNCH.EQ. 1.0R.NPUNCH.EQ.3)  THEN 
WRITE(8,923) 

WRITE(8,924)  ' EPSLON ' , ITRATE, ' : ’ 
DO  601  I=1,NW-1 

WRITE (8, 921)  ERSAVE(I) , ' 

601  CONTINUE 

WRITE (8, 922)  ER3AVE(mj) 
WRITE(8,923) 

END  IF 


C 

WRITE  (6,908)  IDENT 

CALL  PPLOT (PLO(5) , 137, PUP (5) , 137,NW-4, 1, 6HPLO-UP, 137) 
IF(ITRATE.EQ.O)  GO  TO  501 
ITRATE  =  ITRATE  -  1 

CALL  SMOOTH (ERROR, ERRF, WINW, NW,WIND, PLAB) 

GO  TO  59 

501  WRITE(6,908)  IDENT 

CALL  APLOT(PLO(5) , KR, 1 , 2, NW-4 , 1, 2, PLAB ( 5 ) , DUMMY ) 
WRITE(6,915) 

DO  800  J=1,NC 
SUMl  =0.0 
SUM2  =0.0 
DO  700  1=1, NR 

SUMl  =  SUMl  +  (HT(I,J)  *  S(I))  **  2 
SUM2  =  SUM2  +  HT(1,J)  *  A(1,J) 

700  CONTINUE 

ERRJl  =  SORT (SUMl) 

DIFFR  =  1.0  -  SUM2 
UHAQJ  =  ABS (DIFFR)  *  Q(J) 

WRITE ( 6 , 91 6 ) J, ELAB (J),X(J),Q(J), ERRJl , UHAQJ, IDENTY 
I F ( NPUNCH . EQ . 2 . OR . N  PUNCH . EQ . 3 )  THEN 

WRITE (7, 976) J, ELAB (J) ,X(J) ,Q( J) , ERRJl, UHAQJ, IDENTY 
END  IF 
800  CONTINUE 


C 


CLOSE  (2) 

CLOSE  (3) 

CLOSE  (5) 

CLOSE  (6) 

CLOSE  (7) 

CLOSE  (8) 

901  FORMAT (8110) 

902  FORMAT ( 6X, ' STRETCH 
+' NUMBER  OF  BITS 

+' NUMBER  OF  ROWS 
+' NUMBER  OF  COLUMNS 
+' NUMBER  OF  WINDOWS 
+' FIRST  SPECTRUM  NUMBER 


=',F5.3/6X, 

=',I5/6X, 

=',I5/6X, 

=',I5/6X, 

=',I5/6X, 

=’,I5/6X, 
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+’ PUNCHING  OPTION  = ' , 15, 6X, ’ ( 0/1/2/3)  --  (NO  PUNCHING  / 

+  SMOOTHED  OUTPUT  /  UNSKOOTHED  OUTPUT  /  BOTH)'/6X, 

+  'NUI<1BER  OF  ITERATIONS  =',I5/6X, 

+’ PERCENTAGE  ERROR  =', I5/6X, ' FIRST  OVERLAP  BIN  =’,I5) 

903  FORMAT(1X,I10,F10.3,4I10,F10.4) 

904  FORMAT(lX,E13.3,I5,3E13.3,F7.2,2E13.3) 

905  FORMAT (20A4) 

906  FORMAT (1H1,20A4) 

907  FORMAT{6X, ' ENERGY' , 14X, 'PLO',10X, 'PUP’,9X, ’PAVE',6X, 

+'PCT  W',6X, 'ERR1',9X, 'ERR2',/) 

908  FORMAT ( IHl, 'PLO  AND  PUP ' , 6X, 20A4 ) 

909  FORMAT  ( IX,  '****THIS  IS  AN  INVALID  RUN,  ...'.FOR  DIAGNOSTIC  ONLY') 

910  FORMAT(lX,I3,3E12.3,2F12.3,E12.3) 

911  FORMAT  (3X, 'I',6X, 'VLO',10X, 'B',11X, ’S',9X, 'BADJ',9X, 

+'PCTDEV',8X, 'ESD',/) 

912  F0RMAT(1X,I3,4E12.3,2F12.3) 

913  FORMAT  (IHl, '  INDICATOR  NO.  1  =’,F12.3,*  (CONDITION)',//, 

+'  INDICATOR  NO.  2  =',F12.3,'  (CONSISTENT  WITH  NON-NEG. )',// , 

+'  INDICATOR  NO.  3  =',F12.3,'  (RESIDUALS)',//, 

+'  SUM  OF  INDICATORS  =',F12.3. 

+'  (SHOULD  BE  LESS  THAN  5.000  FOR  VALID  RUN  —  SEE  MANUAL) ') 

914  FORMAT(10X,E10.4,I5,5X,E10.4,5X,E10.4,I5,5X,E10.4) 

915  FORMAT ( IHl, 8X, 'J',8X, 'ELAB',9X, 'X',12X, 'Q',9X, 'ERRJ1',6X, 'UHAQJ',6 
+X, 'SPECTRUM  NUMBER' /IHO) 

916  FORMAT(5X,I5,4X,E11.4,2X,E11.4,2X,E11.4,2X,E11.4,2X,E11.4,2X,I3) 

976  FORMAT(2X,I5,2X,E11.4,2X,E11.4,2X,E11.4,2X,E11.4,2X,E11.4,2X,I3) 
918  FORMAT (' -  TAW=',F8.3, '  - ’,/) 

920  FORMAT (IX, A3, ll,Alj 

921  FORMAT(1X,E10.4,A1) 

922  FORMAT(1X,E10.4) 

923  FORMAT (IX) 

924  FORMAT (IX, A6, II, Al) 

C 

STOP 

END 

C  SUBROUTINE  BSREAD 

C  This  routine  assumes  that  B  and  S  are  contiguous  in  storage. 

C  -  IL  -  is  the  first  bin  of  the  8-bin  overlap 

C - 

SUBROUTINE  BSREAD ( B, S , MR, NR, NC, IDENT, LLO, LUP, VLO, VUP, FORTE, IL) 

C 

DIMENSION  B(NR) ,S(NR) ,IDENT(20) 

DIMENSION  VLO (NR) , VUP (NR)  ,  LLO (NR) , LUP (NR) , FORTE (NR) 

C 

CALL  TOWCBN  (B, S, NR, IDENT, LLO, LUP, VLO, VUP, FORTE) 

1  SCOFF=.01 
FUDGE=1 . 0 
DO  137  1=1, NR 

S ( I ) =AMAX1 { S ( I )  ,  SCOFF*B ( I ) ) 

137  CONTINUE 

COVH=B ( IL) +B ( IL+2 ) +B ( IL+4 ) +B ( IL+6) +1 . E-20 
SOVH=SQRT (S (IL) **2+S (IL+2) **2+S (IL+4 ) **2+S (IL+6) **2) 

RERRH=SOVH/COVH 

COVL=B(IL+l) +B(IL+3) +B ( IL+5) +B (IL+7) +1 . E-20 

SOVL=SQRT (S (IL+1) **2+S (IL+3) **2+S (IL+5) **2+S (IL+7) **2) 

RERRL=SOVL/COVL 

FUDGE=COVL/COVK 


RELFU=ABS ( 1 . -FUDGE) 

C 

C  Decide  on  strategy. 

C 

12  IF{RERRH-RELFU) 13, 13,  14 

13  IF(RERRL-RELFU) 16, 16, 14 

14  FUDGE=1.0 
GO  TO  34 

16  FUDGEH=1 . + (COVL-COVH) /COVH*SOVL/ (SOVH+SO’/L) 

17  FUDGEL=1.+(COVH-COVL)/COVL*SOVH/(SOVH+SOVL) 

DO  20  1=1, IL 

B(I)=B(I) *FUDGEH 
20  S(I)=S(I)*FUDGEH 
N1=IL  +  2 
N2=IL  +  6 
DO  30  I=N1,N2,2 
B(I)=B(I) *FUDGEH 
30  S(I)=S{I)*FUDGEH 
N1=IL  +  7 
DO  32  I=N1,NR 
B(I)=B(I) ‘FUDGED 

32  S(I)=S(I)*FUDGEL 
N1=IL  +  1 
M2=IL  +  5 

DO  33  I=N1,N2,2 
B{I)=B(I) ‘FUDGED 

33  S(I)  =  S(I)*FUDGED 

90  FORMAT  (IHl,  45HCOUNTS  AND  STANDARD  DEVIATION  IN  GROUPED  BINS 
+  ,5X,7HSCOFF  =, F8 . 4 , 5X, 7HFUDGE  =, F8 . 4 , 5X, 7HREDFU  =,F8.4) 

34  WRITE  (6, 90) SCOFF, FUDGE, REDFU 

CALD  APDOT(B,MR, 1,2, NR, 1,2, VDO, DUMMY) 

RETURN 

END 

C  SUBROUTINE  WREAD 

C - 

SUBROUTINE  WREAD (W, NO, NW, NTRY, EDAB, PDAB, PCTWID) 

C 

DIMENSION  W(NC) , PDAB (NW) , PCTWID (NW) , EDAB (NC) 

COMMON  /STORE/  SIGJ(llO) 

C 

GAUSS(X)  =  0.398942280414  *  EXP(-0.5  *  X**2) 

IF(NTRY.EQ. 1)  CADD  WCADC (NC, NW, EDAB, PDAB, PCTWID) 

DO  10  J=1,NC 

DIFOS  =*  ABS  (EDAB  (J) -PDAB  (NTRY)  )  /  SIGJ(J) 

IF(DIFOS  -  5.5)  3,2,2 

2  W(J)  =  0.0 
GO  TO  10 

3  W(J)  =  l.OE+6  *  GAUSS(DIFOS)  /  SIGJ{J) 

10  CONTINUE 

RETURN 

END 

C  SUBROUTINE  WCADC 


SUBROUTINE  WCADC (NC,NW, EDAB, PDAB, WINW) 
DIMENSION  EDAB(NC) ,PDAB(NW) , WINW (NW) , I ( 120) 
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COMMON  /STORE/  SIGJ(llO) 

C 

DO  10  J=1,NC 
KJ)  =  0 
DO  2  K=1,NW 

IF(ELAB(J)  -  PLA3(K))  1,2,2 

1  I(J)  =  K 
GO  TO  10 

2  CONTINUE 
10  CONTINUE 

DO  20  J=1,NC 
II  =  I(J) 

IF(II.EQ.l)  GO  TO  30 
IF(II.EQ.O)  GO  TO  40 

W  =  WINW(II)  -  (WINW(II)-WINW(II-l) )* (PLAB(II)-EIAB(J) )/(PIAB(II)- 
+  PLAB(Il-l)) 

GO  TO  20 

30  W  •=  WINW(l)-(V;iNW(2)-WINW(l)  )*{PLAB(1)-ELAB(J))/(PLAB(2)-PLAB(1)  ) 
GO  TO  20 

40  W  =  WINW(NW)  +  (WINW(NW) -WINW(NW~1) ) * (ELAB(J)  -PLAB(NW))  / 

+  (FUAB(NW) -PLAB(NW-l) ) 

20  SIGJ(J)  =  (0.01  *  W  *  ELAB(J))  /  2.35482004502 
RETURN 
END 

C  SUBROUTINE  TOWCBN 

C - 

SUBROUTINE  TOWCBN  (GC,  S,  NR, IDENT, LLO, LUP,VLO,VUP, FORTE) 

C 

CHARACTER*6  TiTLEl, TITLE2 

DIMENSION  C(5000) ,GC(NR) ,S (NR) ,LLO(NR) ,LUP(NR) 

DIMENSION  VXO(NR) ,VUP(NR) ,FORTE(NR) 

DIMENSION  CRE(5000) ,BRE(5000) 

DIMENSION  B(5000) ,UC(4) ,UB(4) 

DIMENSION  FGAIN ( 4 ) , BGAIN ( 4 ) , GNOM ( 4 ) 

DIMENSION  FAF(19) 

DIMENSION  TITLE1(4),TITLE2 (4) 

DIMENSION  FBCORR(20) , IDENT(20) 

EQUIVALENCE  (FBCORR ( 1 ) , UC ( 1 ) ) , (FBCORR(5) ,UB(1) ) , 

+(FBCORR(9) ,FGAIN(1) ) , ( FBCORR ( 13) , BGAIN ( 1 ) ) , ( FBCORR ( 17 ), GNOM ( 1 ) ) 

C 

READ  (3,909)  FAF 
909  FORMAT (19F4. 3) 

READ  (3,900)  WPOINT,NGRPS,NGAIN,NPUNCH 

900  FORMAT (416) 

READ (5, 901) (UC(I),1=1,4), (UB (I ) , 1=1, 4 ) 

+  , (FGAIN (I) ,1=1,4) , (BGAIN (I), 1=1, 4) 

READ(3,901) (GNOM(I) ,1=1,4) 

901  FORMAT (IX, 4F10.0) 

CALL  BlNRD(NGRPS,LLO,LUP,NBINS) 

WRITE  (6,400) 

400  FORMATdHl,'  RUN  F  LIVE  TIMES  B  LIVE  TIMES  ’, 

+'F  GAIN  B  GAIN  NOM  GAIN’) 

CALL  FPRINT  (UC,  4 , 4 , 5, 3, 0, 1 , 0) 

CALL  CREAD  (C, 5000, TITLEl , TITLE2 ,NERR,NPOINT, IDENT) 

DO  30  K  =  1,NGAIN 
KK  =  1  +  K*1000 
WRITE  (6,904)  TiTLEl (K),K 
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30  CALL  CPRINT  (C ( XK) , NPOINT) 

CALL  CREAD  ( B, 5000, TITLEl , TITLE2 , NERR, NPOINT/ IDENT) 

DO  31  K=1,NGAIN 
KK  =  1 ^K*1000 

26  WRITE  (6,905)  TITLE1(K),K 

31  CALL  CPRINT (B(KK) , NPOINT) 

DO  320  K=1,NGAIN 

KK  =  1  +  K*1000 

DO  321  1=1, NPOINT 
VVLO=FLOAT(I-l) *GNOM(K) 

VVUP=FLOAT(I) *GNOM(K) 

J  =  KK+1-1 

CRE(J)=REBIN(C(KK)  ,  1,  NPOINT,  0,  FGAIN  (K)  ,  1 . 0,  WLO,  WUP) 
BRE(J)=REBIN  (B(KK)  ,  1, NPOINT,  0,BGAIN(K)  ,  1 . 0,  WLO,  WUP) 

321  CONTINUE 
320  CONTINUE 

DO  322  K=1,NGAIN 

KK  =  1+K*1000 

WRITE  (6,908)  TITLE1(K),K 

322  CALL  CPRINT (CRE(KK) /NPOINT) 

DO  323  K=1,NGAIN 
KK=1+K*1000 

325  WRITE  (6/910)  TITLE1(K),K 

323  CALL  CPRINT (BRE(KK) /NPOINT) 

DO  1999  I=1/NBINS 

GC ( I ) =0 . 0 
VAR  =0.0 
Ll=LLO(I) 

L2=LUP(I) 

LL  =  Ll/1000 

LLLO  =  LI  -  LL*1000 

LLUP  =  L2  -  LL*1000 

VLO(I)  =  FLOAT ( LLLO- 1 ) *GNOM(LL) 

VUP(I)  =  FLOAT (LLUP) *GNOM(LL) 

BINW  =  VUP(I)  -  VLO(I) 

DO  1990  L=L1,L2 

GC (I ) =GC (I ) +CRE (L) /UC (LL) -BRE (L) /UB (LL) 

1998  VAR=VAR+ (CRE (L) +. 1) /UC (LL) **2+ (BRE (L) + . 1 ) /UB ( LL) **2 
GC(I)  =  GC(I) /BINW 

1999  S(I)  =  SQRT (VAR) /BINW 
DO  2  I=1/NBINS 

V=0.5* (VLO(I)+VUP(I) ) 

IF(V-O.IO)  4/8/8 

4  IF(V-O.Ol)  5/  6 

5  FORTE(I)  =  FAF(l) 

GC(I)  =0.0 

S(I)  =  0.0 
GO  TO  2 

8  FORTE (I)  =1.0 

GO  TO  7 

6  Z  =  200.0* (V-0.01) 

J=Z 

Z  =  Z-FLOAT(J) 

FORTE(I)=(1.0-Z)*FAF(J+l)+Z*FAF(J+2) 

7  CONTINUE 

GC(I)  =  GC (I) /FORTE (I) 

S(I)  =  S(I)/FORTE(I)  +  0.2* (l.O-FORTE(I) ) *GC(I) 

2  CONTINUE 
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WRITE  (6,906)  TITLEl(l) 

906  FORI-IAT  (IHl,  lOX, A6,  5X,  19HGROUPED  COUNTS/UNIT, /,  IH  ) 

WRITE  (6,401) 

401  FORMAT  (IH  , 6X, 4HBINS, 6X, 6HCOUNTS, IIX, 5HERROR, 6X, 6HLOCHAN, 2X, 
+6HUPCHAN,6X,7HLO  P.H.,9X,7HUP  P. H.  ,  7X, 6HF. A. F. ) 

DO  50  I=1,NBINS 

50  WRITE  (6,934)  I , GC ( I ) , S ( I ) , LLO ( I ) , LUP ( I ) , VLO ( I ) , VUP ( I ) , FORTE ( I ) 

934  FORMAT  (IH  , 18, 2E16 . 5, 218, 2E16 . 5, FIO . 5) 

999  RETURN 

904  FORMAT  ( IHl , lOX, A6, 5X, 2 8HRAW  COUNTS  FOREGROUND  GAIN  =,I2,/,1H  ) 

905  FORMAT  ( IHl , lOX, A6, 5X, 28HRAW  COUNTS  BACKGROUND  GAIN  =, 12 , / , IH  ) 

906  FORMAT ( IHI, 10X,A6,  5X,  30HREBIN  COUNTS  FOREGROUND  GAIN  =,I2,/,1K  ) 
910  FORMAT ( IHl, 10X,A6,5X,30HREBIN  COUNTS  BACKGROUND  GAIN  =,I2,/,1H  ) 

END 

C  SUBROUTINE  MATIN 

C  This  subroutine  reads  data  into  an  array. 

C  A  format  card  must  precede  the  data,  i.e.  FORMAT  I  NN  (.,.),  with: 

C  'FORMAT'  in  col.  1-6, 

C  I  =  'ROW'  for  row  format,  I  =  'COL'  for  column  format, 

C  NN  =  number  of  data  words  per  card, 

C  (...)=  format  specification. 

C  Each  data  card  must  have  two  index  numbers  indicating  the  matrix 

C  position  of  the  first  data  word.  The  reading  of  cards  is 

C  completed  by  one  blan){  card.  This  routine  has  a  provision  for 

C  changing  formats  in  midstream.  A  card  with  both  indices  negative 

C  implies  that  a  format  card  follows. 

c - 

SUBROUTINE  MATIN (A, MR, NR, NC/ 

C 

DIMENSION  A(MR,NC),FMT(17),XLIST(78) 

C 

CALL  EQUIV(ICHECK,4HFORM) 

CALL  EQUIV(IROW,4HROW  ) 

NERROR=0 
DO  4  I  =  1, MR 
DO  4  J  =  1,NC 

4  A(I,J)  =  0.0 
C  ***Read  format  card 

5  READ  (3,900)  IWORD, IROC, NN, FMT 
IF  (IWORD/2-ICHECK/2)  8,10,8 

8  WRITE  (6,904) 

GO  TO  888 

C  ***Read  data  card 

10  READ  (3,  FMT)  I,  J,  (XLIST  (K)  ,  K=1,NN) 

C  ***Chec)c  for  new  format  card 

32  IF(I)50,11,11 
50  IF(J)5,721,701 

C  ***Checl<  for  blan)c  card  to  end  list 

11  IF  (I  +J)  12,888,12 

12  1  =  MAX0(I,1) 

J  =  MAX0(J,1) 

IF  (IROC-IROW)  213,13,213 
C  ***  Row  format  *** 

13  IF  (I-NR)  14,14,701 

701  IF(NERROR-50)  702,702,10 

702  WRITE  (6,901)  I, J,MR,NR,NC, (XLIST (Kl) , Kl=l, 4 ) 

NERROR=NERROR+l 
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o  o  non 


GO  TO  10 

14  DO  20  JJ  =  1,NN 

JJJ  =  J  +  JJ  -  1 
IF  (XLIST(JJ) ) 16,20, 16 

16  IF  (JJJ  -  NO  17,17,703 

703  IF(NERROR-50)  704,704,10 

704  WRITE  (6,902)  I, JJJ,MR, NR, NC, (XLIST (Kl) , Kl=l, 4 ) 

NERROR=NERROR+l 

GO  TO  10 

17  IF  (A(I,JJJ))  705,18,705 

705  IF(NERROR-50)  706,706,20 

706  WRITE  (6,903)  I , JJJ, MR, NR, NC, (XLIST (Kl) , Kl=l, 4 ) 

NERROR=NERROR+l 

GO  TO  20 

18  A(I,JJJ)  =  XLIST(JJ) 

20  CONTINUE 

21  GO  TO  10 

C  *■*■*  Column  format  *** 

213  IF  (J-NC)  214,214,721 

721  IF(NERROR-50)  722,722,10 

722  WRITE  (6,  902)  I , J, MR,  NR, NC,  (XLIST (Kl )  ,  Kl=l,  4 ) 

NERROR=NERROR+l 

GO  TO  10 

214  DO  220  II  =  1,NN 

III  =  I  +  II  -  1 
IF  (XLIST(II))  216,220,216 

216  IF  (.TII-NR)  217,217,723 

723  IF(NERROR-50)  724,724,10 

724  WRITE  (6,901)  III , J,MR, NR, NC, (XLIST (Kl ) , Kl=l, 4 ) 

NERROR=NERROR+l 

GO  TO  10 

217  IF  (A(III,J))  725,218,725 

725  IF(NERROR-50)  726,726,220 

726  WRITE  (6,903)  III, J,MR, NR, NC, (XLIST (Kl ) , Kl=l, 4 ) 

NERROR=NERROR+l 

GO  TO  220 

218  A(II1,J)  XLIST(II) 

220  CONTINUE 

GO  TO  10 
888  RETURN 

900  F0RMAT(1A4,3X,1A3,1X,I2,1X,16A4,1A2) 

901  FORMAT  (22H0ILLEGAL  1ST  INDEX,  I=, 13, 2X, 2HJ=, 13, 2X, 3HMR=, 13, 2X, 

+  3HNi<=,  I3,2X,  3HNC=,  13,  5H  **,  4E12 . 3,  3H.  .  . ) 

902  FORMPT  (22H0ILLEGAL  2ND  INDEX,  I=, 13, 2X, 2HJ=, 13, 2X, 3HMR=, 1 3, 2X, 

+  3HNR=,I3,2X,3HNC=,I3,5H  *♦, 4E12 . 3, 3H. . . ) 

903  FORMAT  (22HODUPLICATION  ERROR,  I=, 13, 2X, 2HJ=, 13, 2X, 3HMR=, 13, 2X, 

+  3HNR=,I3,2X, 3HNC=,I3,5H  **, 4E12 . 3, 3H. . . ) 

904  FORMAT (1H0,51H***  PROGRAM  EXPECTS  FORMAT  CARD  AT  THIS  POINT  ***) 
END 


SUBROUTINE  EQUIV 


SUBROUTINE  EQUIV (X,Y)' 

X  =  Y 

RETURN 

END 


SUBROUTINE  AREAD 
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n  o 


c 


SUBROUTINE  AREAD {A, MR, NR, NC, NW) 

DIMENSION  A(l) 

CALL  MATIN (A, MR, NR, NC) 

RETURN 

END 

C  SUBROUTINE  BINRD 

C  Reads  bin  card(s)  with  (KLO (I) , KUP (I ) , KDEL ( I ) ) , 1=1 , NGRRS 
C  and  converts  to  expanded  bin  list  (LLO ( J) , LUP ( J) ) , J=l, NBINS . 

C  Needs  NGRPS  from  calling  program. 

C  Returns  LLO, LUP, NBINS. 

C - 

SUBROUTINE  BINRD (NGRPS, LLO, LUP, NBINS) 

C 

DIMENSION  LLO (7777) , LUP (7777) , KUP (40) , KLO (40) , KDEL ( 40) 

Read  binning  card(s),  10  max 
READ (3, 901)  (KLO(I) ,KUP(I) ,KDEL(I),I=1,NGRPS) 

C  Write  out  binning  table 

WRITE  (6,902) 

DO  3  M=l, NGRPS 
AUP  =  KUP(M) 

ALO  =  KLO(M) 

ADEL  =  KDEL(M) 

BIN  =  (AUP-ALO+1.0)/ADEL 
WRITE  (6,  903)  M,  KLO (M)  ,  KUP (M) , KDEL (M) , BIN 
3  CONTINUE 
C  Convert  table 

NBINS=0 

DO  20  1=1, NGRPS 
Kl=KLO(I) 

K2=KUP(I) 

K3=KDEL(I) 

K3  =  MAX0(K3, 1) 

DO  10  K=K1,K2,K3 
NBINS=NBINS+1 
LLO (NBINS) =K 
LUP (NBINS) =K+K3-1 
10  CONTINUE 
20  CONTINUE 

C  Write  out  computed  number  of  bins  I 

WRITE  (6,904)  NBINS  1 

RETURN  i 

901  FORMAT (1216)  1 

902  FORMAT ( IHl , 13HBINNING  TABLE/6HOGROUP, 4X, 6HtOEDGE, 4X, 6HUPEDGE, 4X, 8H 

+CHAN/BIN,5X,4HBINS)  i 

903  FORMATdH  ,  13,  6X,  15,  5X,  15,  8X,  13,  7X,  F5. 1)  \ 

904  FORMAT  (42X,4H - /25X,  14HTOTAL  NO.  BINS,4X,!l3) 

END  1 

C  SUBROUTINE  FPRINT  1 

C  Output  routine  for  1  or  2  dimensional  arrays!. 

C  I 

C  K  =  1,  FORMAT  112 

C  K  =  2,  FORMAT  F12 . 4 

C  K  =  3,  FORMAT  E12.4 

C  K  =  4,  FORMAT  G12 . 4 
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C  MR  =  first  dimenstion  number  of  array 

C  L  =  0  for  no  column  headings 

C  N  =  0  for  return  without  action 

C  LL  =  1  for  placing  fixed  point  numbers  in  G  format 

C 

C  Compiled  in  FORTRAN  360  for  system  360  on  25  July  1966. 

C - 

SUBROUTINC  FPRINT  (A, MR, NR, NC, K, L, N, LL) 

C 

DIMENSION  A  {MR,  888)  ,  FMT’(5),  FMT3  ( 9)  ,  FMT4  ( 4  ) 

DATA  FMT/4H(I6, ,4H2X,  ,  4H  10, 4H  ,4H4)  / 

DATA  FMT3/4H  1,4H  2,4H  3, 4H  4,4H  5, 4H  6, 4H  7, 

+  4H  8,4H  9/ 

DATA  FMT4/4HI12) ,4HF12.,4HE12.,4HG12./ 

C 

4  IF  (N)  1,99,1 
1  NTOU  =  6 
C***  Set  format 

DATA  FMT(3)/4H  10/ 

FMT(4)  =  FMT4(K) 

C 

JUMP  =  0 

18  WRITE  (NTOU, 902) 

19  NFULLB  =NC/10 

20  IF  (NFULLB)  21,40,21 

21  NSKIP  =  MAXO (1, 60/ (NR+4) ) 

DO  39  M=l, NFULLB 

NCI  =  10*M  -  9 

NC2  =  10*M 

GO  TO  80 

37  IF  (M  -  NSKIP* (M/NSKIP) )  39,38,39 

38  WRITE  (NTOU, 903) 

39  CONTINUE 

40  NCLEFT  =  NC  -  10*NrJLLB 

41  IF  (NCLEFT)  42,99,42 

42  NCI  =  10*NFULLB  +  1 

NC2  10*NFULLB  +  NCLEFT 
FMT(3)  =  FMT3 (NCLEFT) 

JUMP  =  1 

80  WRITE  (NTOU, 902) 

IF  (L)  81,83,81 

81  WRITE  (NTOU, 901) ( (J),J=NC1,NC2) 

82  WRITE  (NTOU, 902) 

83  WRITE  (NTOU,FMT) (I, (A(I,J),J=NC1,NC2),I=1,NR) 

IF  (JUMP)  99,37,99 

99  RETURN 

901  FORMAT (IH  ,3X, 10112) 

902  FORMAT (IX) 

903  FORMAT (/IHl) 

END 

C==============================:===================================; 

C  SUBROUTINE  CREAD 

C  Reads  spectral  data  from  input  file  SPECTRA. IN. 

C  Format  lines  in  SPECTRA. IN  must  have  the  form: 

C  col  1-6  DYNMAT 

C  col  14-31  10  (1A6,I4,10F7.0) 

C - 

SUBROUTINE  CREAD  (C, N, TITLEl , TITLE2, NERR, NPOINT, IDENT) 
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c 

CHARACTER* 6  TITLEl , TITLE2 

DIMENSION  C(777) ,CARD(20) ,FMT(20) ,IDENT(20) 
DIMENSIC;:  TITLEl  (4)  ,TITLE2  (4) 

C 

NERR  =  0 
KK  =  0 

DO  30  K  =  1,N 
C{K)  =  0.0 
30  CONTINUE 

1  READ  (5,901)  CARD(l) ,IREL,NWORD,  {EMT(J) , J=l,10) 
NCHAN  =0 

KK  =  KK+1 

CALL  EQUIV  (X,4HDYNM) 

IF  (CARD(l)/2.0-X/2.0)  3,2,3 
3  CALL  EQUIV(X,4H  ) 

IF  (CARD(l)/2.0-X/2.0)  2,88,2 

2  IF  (IREL  .EQ.  1000)  THEN 

WRITE  (6,903/  CARD ( 1 ), IREL, NWORD,  (EW  (J) , J=l, 10) 
ELSE 

WRITE(6,906)  CARD ( 1 ), IREL, NWORD, (FMT ( J) , J=l, 10) 
END  IF 

READ  (5,900)  ( IDENT ( 0 ) , J=1 , 20 ) 

WRITE  (6,904)  ( IDEN f ( J) , J=l, 20! 

10  READ(5,FKr)  TITLEl ( KK) , 10, (CARD ( J! , J=l, NWORD) 
WRITE(6,905)  TITLEl (KK) , 10, (CARD (J) , J=l, NWORD) 
TITLE2 (KK) =TITLE1 (KK) 

10  =  lO+IREL 
DO  17  I  =  1, NWORD 
NCHAN  =  NCHAN+1 
C 

C  * 

C  . INDEX  OK . 

C  NO  YES 


C  ... .WORD  ZERO. 


C 

NO 

YES 

STORE 

C 

* 

* 

C 

ALARM 

* 

★ 

C 

* 

* 

★ 

C 

. CONTINUE . 

c 

. CONTINUE . 

C 

■k 

C 

II  = 

10  +  I  -  1 

12 

IF 

(II) 

14,14,13 

13 

IF 

(II-N) 

16, 16, 14 

14 

IF 

(CARD (I) ) 

99, 17, 99 

16 

C(II) 

=  CARD (I) 

17 

CONTINUE 

IF  (NCHAN-NPOINT)  10,1,10 
88  RETURN 
99  NERR  =  NERR+1 


WRITE  (6,902)  N,CARD 
GO  TO  10 

900  FORMAT  (1X,20A4) 

901  FORMAT  ( IX, A6, IX, 16, 13, IX, 10A4 ) 

902  FORMAT( 'O’, ’***  CREAD  ERROR,  N  =  ',16,'  CARD  IS  **’,20  A4,  ’**’) 

903  FORMAT  ( IHl, A6, 2X, 16, 13, IX, 10A4 ) 
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904  FORMAT  (IH  ,2CA4) 

905  FORMAT  (IH  , 1A6, 1 1 6, 1 0F8 . 0, 1A6) 

906  B'ORMAT  ( 1H0,A6,  2X,  16,  13,  IX,  10A4  ) 
END 


SUBROUTINE  CPRINT 


SUBROUTINE  CPRINT (C,NP) 

C 

DIMENSION  C(7777) ,XBUE(8) ,IBUE{8) 

C 

TCOUNT  =0.0 
NPAGES  =  1  +  (Nr-l)/400 
DO  10  K=l, NPAGES 
IF  (K-1)  5,5,4 

4  WRITE  (6,904) 

5  NPOINT  =  MINO (400,NP- (K-1) *400) 

NROWS  =  MINO (NPOINT,  50) 

DO  7  1=1, NROWS 

NCOLS=(NPOINT/50) +MIN0 '.1,  (NPOINT-50*  (NPOINT/50)  ) /I ) 

NCOLS  =  MAXO ( l,NCOL3) 

DO  8  J=l, NCOLS 

NC=I+ ( J-1) *50+ (K-1) *400 
XBUE(J)=C(NC) 

TCOUNT  =  TCOUNT+XBUE ( J) 

IBUE(J)=NC 
8  CONTINUE 

WRITE  (6,  905)  ( IBUE (L) , XBUE (L) , L=l,  NCOLS) 

7  CONTINUE 
]0  CONTINUE 

WRITE  (6,906)  TCOUNT 
RETURN 

904  FORMAT (IHl) 

905  FORMATdH  ,  8  (15,  F9. 1, 2X)  ) 

906  FORMAT  ( IHO/lHO, 24HTOTAL  NUMBER  OF  COUNTS  =,F12.1) 

END 

C  FUNCTION  REBIN 

C  The  vector  A  gives  the  components  fo  a  histogram.  The  bottom 
C  edge  of  the  lower  bar  is  a  RiiOLO,  and  the  width  of  the  bars  is 

C  RHODEL.  This  routine  adds  all  teh  bins  between  VLO  and  VUP  and 

C  splits  any  fractional  bins  linearly.  For  greater  flexibility, 

C  V  =  C*RHO,  where  C  is  an  input  parameter. 

C-  - - 

FUNCTION  REBIN (A, MR, N, RHOLO, RHODEL, C, VLO, VUP ) 

C 

DIMENSION  A(777) 

C 

1  RLO=AMAXl (VLO/ C,  RHOLO) 

2  RUP=AMAX1 (VUP/C, RHOLO) 

C 

10  XI  =  2.0  +  (RLO- RHOLO) /RHODEL 
II  =  XI 

II  =  MIN0(I1,N+1) 

XI  =  Xl-FLOAT(Il) 

XI  =  AMINl (XI, 1.0) 

C 

15  X2  =  1.0  +  (RUP- RHOLO) /RHODEL 
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12  =  X2 

12  =  MIN0(I2,N) 

X2  =  X2-FLOAT(I2) 

X2  =  AMINl (X2, 1. 0) 

C 

20  LI  =  ( (Il-l)  -  1)*MR  +  1 
L2  =  (12-1) *MR  +  1 

REBIN  =  (1.0-X1)*A(L1)  -  ( 1 . 0-X2) *A(L2) 

C 

IF  (I2-I1+1)  40,40,30 

30  DO  37  I  =  11,12 

L3  =  (I-1)*MR  +1 
REBIN  =  REBIN  +  A{L3) 

37  CONTINUE 
C 

40  RETURN 
END 

C  SUBROUTINE  GINV 

C  Upon  entry,  A  =  a  matrix  with  NR  rows  and  NC  columns. 

C  After  the  return,  the  original  A  is  destroyed  and  the  array  A 
C  contains  the  transpose  of  T,  the  transformation  matrix 
C  (X  =  T*B)  which  solves  the  pari  of  equations: 

C  A  *  X  =  B 

C  and 

C  TAU  I  *  X  =  0 

C  in  the  least  squares  sense. 

C  MR=  1st  dimension  no.  of  array  A  in  the  calling  program. 

C  TAU  is  a  nonnegative  constant  which  controls  the  error 
C  propagation  of  the  transformation  matrix.  If  TAU=0,  then  the 
C  resulting  transformation  matrix  si  the  ordinary  least  squares 
C  transformation  matrix  (the  inverse  if  NR=NC) .  ATEMP  and 

C  U  are  used  for  working  space  by  the  algorithm  and  do  not 

C  necessarily  contain  any  relevant  numbers  at  the  conclusion. 

C  ATEMP  must  be  dimensioned  at  least  NC  by  the  calling  program. 

C  U  must  be  dimensioned  at  least  NC*NC  by  the  calling  program. 

C  TAU  is  a  nonnegative  constant  which  controls  the  error. 

C  Number  of  multiplications  =  NC**2  (5/2  +  2/3  NC) . 

C - 

SUBROUTINE  GINV (A, MR, NR, NC, TAU, U, ATEMP) 

C 

DIMENSION  A(MR,r:C)  ,U(NC,NC)  ,  ATEMP  (NC) 

C 

TAUSQ  =  TAU  **2 
C  Place  unit  m.atrix  in  U 
DO  5  1=1, NC 

DO  4  J  =1,NC 

4  U(I,J)  =  0.0 

5  U(I,I)  =  1.0 
C 

C  Orthogonalize  combined  matrix  (A  above  U)  by  Gramm-Schmidt-Hilbert 
C  method  with  first  NR  rows  weighted  with  1  and  the  other  NC  rows 
C  weighted  with  1/TAU.  Then  reorthogonalize  to  lessen 

C  roundoff  error. 

C 

DO  20  1=1, NC 
II  =  I-l 
IF  (II)  2,11,2 
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DO  10  LL=1,2 
DO  10  J=1,II 

DOT  =  CDOT(A(l,I) ,A(1,J),NR) 

DOT2  =  CDOT(U(l,I),U(l,J),J) 

DOT  =  DOT  +  DOT2  *  TAUSQ 
CALL  CLIM  (U(1,I),U(1,J),DOT,J,0) 

10  CALL  CLIM  (A(I,I),A(1,J),DOT,NR,0) 

Normalize  the  column  *1*  of  the  combined  matrix 


11  DOT  =0.0 
DOT2  =0.0 

DO  12  K  =  1,NR 

12  DOT  =  A(K,I)  *  A(K,I)  +  DOT 

DO  14  K=1,I 

14  DOT2  =U(K,I)  *U<K,1)  +  DOT2 

DOT  =  DOT  +  DOT2  *  TAUSQ 
DOT  =  SQRT(DOT) 

DO  17  K=1,I 

17  U(K,I)  =  U(K,I)/DOT 

DO  19  K=1,NR 

19  A(K, I)  =  A (K, I) /DOT 

20  CONTINUE 

Calculation  of  the  transpose  of  the  transformation 
matrix  T(trans.)  =  A  +  U(trans.) 


DO  50  1=1, NR 

DO  45  J=1,NC 

ATEMP(J)  =  RDOT(A(I, J) ,MR,U(J,J),NC,NC-J+1) 
45  CONTINUE 

DO  50  J=1,NC 

A(I,J)  =■  ATEMP(J) 

50  CONTINUE 
RETURN 
END 
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c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 
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SUBROUTINE  APLOT 

This  version  plots  3  plus  and  2  minus  cycles  on  a  p-log  scale. 

row  plots:  NDC  =  1 

NDR  “  first  dimension  of  A 
NCR  =  number  of  rows  in  plot 
NRC  =  number  of  columns  in  plot 

1 

NDC  =  first  dimension  of  A 
NDR  =  1 

NCR  =  number  of  columns  in  plot 
NRC  =  number  of  rows  in  plot 

Every  Nth  tow  or  column  of  matrix  A  is  plotted.  NPPP  specifies 
the  number 'i  of  plots  per  page  (5  max).  X  is  the  absicsa 
labelling  array.  NT  is  the  number  of  the  output  tape. 


column  plots; 


SUBROUTINE  ij^PLOT  (A,  NDC,  NDR, NCR,  NRC, NTH, NPPP,  X, NT) 

DIMENSION  Ai777) ,X(333) ,PLOT(101),ESCALE(11),SAVE(5)  ,SYM{5) 
+,L(5) ,LOC(5) 


ns 


BAKF(Y,F)  =  0,4342944819*ALOG(Y+F) 

FORF(X,F)  =  10.0**X-F 
IF(KKK-14232) 10,30,  10 
10  KKK=14232 

CALL  EQU1V(3YM<1),1H*) 

CALL  EQUIV(SYM(2)  ,  1H$) 

CALL  EQUIV(SYM(3) , 1H=) 

CALL  EQUIV(SYM(4) , 1H+) 

CALL  EQUIV(SYM(5) , IH*) 

DO  15  1=1,101 
15  CALL  EQUIV(PLOT(I)  ,  IH  ) 

DO  20  1=1,101,10 
20  CALL  EQUIV(PLOT(I)  ,  IHl) 

30  INDEX=l 
Ll  =  l 

NPLOT=l+ (NCR-1) /NTH 
40  IF(NPLOT-NPPP)42,44,44 
42  JP=NPLOT 
GO  TO  46 
44  JP=NPPP 
46  NPLOT=NPLOT-NPPP 
JP1=JP+1 
AMAX=1. OE-30 
DO  60  J=1,JP 

L(J)=L1+(J-1) *NTH 
DO  60  1=1, NRC 

IJ=1NDEX+NDR* (I-l) +NDC* { J-1) *NTH 
IF  (AMAXl (A(IJ) ,-10.0*A(IJ) )  -  AMAX)  60,60,50 
50  AMAX  =  AMAXl (A(IJ) ,-10. 0*A(IJ) ) 

60  CONTINUE 

L1=L(NPPP)+NTH 
RSP  =  0 

NSCALE  =  BAKF(AMAX,RSP)+40.95 

NSCALE=NSCALE-40 

FNSM3  =  NSCALE-3 

FO  =  FOKF(FNSM3,RSP) 

FIX=-20 ‘NSCALE 
FIX=FIX+101.5 
FIXN=83. 0-FIX 
SCALE  =  NSCALE 
DO  70  1=1,11 
DIF  =  1-5 

ESCALE(I)  =  FORF(SCALE-3.+.5*ABS(DIF  ),F0) 

70  ESCALE(I)=SIGN(ESCALE(I),DIF) 

IF  (NDR-1)  90,90,100 
90  WRITE  (6,905)  (SYM ( J) , L ( J) , J=l, JP) 

WRITE  (6,902)  ESCALE 
GO  TO  110 

100  WRITE  (6,906)  (SYM(J) ,L(J) , J=l, JP) 

WRITE  (6,903)  ESCALE 

110  DO  130  1=1, NRC 

DO  121  J=1,JP 

IJ=INDEX+NDR* (I-l) +NDC* ( J-1) *NTH 
IF(A(IJ))112,114,114 

112  IJ  =  -20.0*BAKF(F0,-A(1J) )+FlXN 

GO  TO  116 

114  IJ  =  20.0*BAKF(A(IJ) ,F0)+FIX 

116  IF(IJ) 117,117, 118 


116 


117 

IJ=1 

118 

IF!IJ-101) 120, 120, 119 

119 

IJ=101 

120 

10C!J)=IJ 

SAVE!J)=PLOT!IJ) 

121 

PLOT!IJ)=SYM!J) 

WRITE  (6,904)  I,X(I),PLOT 

DO  130  J=1,JP 

JJ=JP1-J 

IJ=LOC(JJ) 

130  PLOT{IJ)=SAVE{JJ) 

INDEX=INDEX+NPPP*NTH*NDC 
IF(NPLOT) 140,  140,  135 
135  WRITE  (6,907) 

GO  TO  40 
140  RETURN 

902  FORMAT! '0',2X,  ' ROW ’ ,  8X,  5E1C . 1 , 3E10 . 1, 3E9 . 1 ) 

903  FORMAT ( '0',2X, 'COL' ,  8X, 5E10. 1, 3E10. 1, 3E9. 1) 

904  FORMAT!'  ' , 13, 2X, ElO. 3, 3X, lOlAl ) 

905  FORMAT!'0',13X,5!12X,A1, '  COL', 13)) 

906  FORMAT!'0’,13X,5!12X,A1,  '  ROW, 13)) 

907  FORMAT  (IHH 
END 

C  SUBROUTINE  PPLOT 

C  Line-by-line  pplot. 

C - - 

SUBROUTINE  PPLOT ! VI ,  MRl ,  V2 , MR2 , N, NO, LABEL, NOPT ) 

DIMENSION  VI (MRl) ,V2 (MR2) ,BUF!102) 

C 

C  Calculate  minim’irn  and  maximum  value  to  be  plotted 
C 

VMAX=-l,0E+35 
VM1N=  l.OE+35 
SUM=0.0 
DO  30  1=1, N 

IF!V1!I)-VMAX) 16,16,14 
14  VMAX=V1!I) 

16  IF!V1 (I) -VMIN) 18,20,20 

18  VMIN=V1!I) 

20  IF!V2 !I)-VMAX)24,24,22 

22  VMAX=V2(I) 

24  IF!V2 (I)-VMIN)26,28,28 

26  VMIN=V2!I) 

28  SUM=SUM+ABS!V1!I)+V2!I)  ) 

30  CONTINUE 
C 

C  Calculate  scale  factor 

C 

TWON=2*N 

AVRGE=SUM/TWON 

NEXP=0 

IF!AVRGE)32,35,32 
32  NEXP=ALOG!AVRGE) *0.4342944819 
35  FEXP=10.0**NEXP 
C 

C  Set  graph  limits 

C 
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GRIDLO=VMIN/FEXP 
GRIDUP=VMAX/FEXP 
DELTA= ( GRI DUP- GRI DLO ) / 1 0  0 . 0 
C 

C  Section  to  place  zero  on  grid  line 

RANGE= 12 0.0* DELTA 
IRANGE^RANGE 
REM= RANGE- FLOAT ( I RANGE ) 
IF(REM)40,70,40 
40  IF(REM-0.5)50,70,60 
50  REM=0.5 
GO  TO  70 
60  REM-1.0 

70  DELTA=(FLOAT(IRANGE)+REM)/100.0 
RATIO=GRIDLO/ (DELTA* 10.0) 
NRATIO=RATIO 
RRATIO=NRATIO 
RE.M-RATIO-RRATIO 
IF(REM)75,80,80 
75  RRATIO=RRATIO-l.  0 
80  GRI DLO=RRATIO* DELTA* 10.0 
GRIDUP=100.0*DELTA+GRIDLO 

Graph  is  now  printed,  line  by  line 

WRITE  (6,900)  LABEL, NO, FEXP 
DO  200  M=1,N 

CALL  EQUIV(BUF(1) ,M) 

MM1=M-1 

IF(M-N)85,140,85 
85  IF(M0D(MM1, 10) ) 90, 140, 90 

90  IF(MOD(MM1,10)-1)180,100,180 

100  DO  110  Ll=l,101 

CALL  EQU1V(BUF(L1+1) , IH  ) 

110  CONTINUE 

DO  120  L2-l,101,10 

CALL  EQUIV ( BUF ( L2+1 ) , IHl ) 

120  CONTINUE 

GO  TO  180 

140  DO  150  L3=l,101 

CALL  EQUIV ( BUF ( L3+1 ), 1H-) 

150  CONTINUE 

DO  160  L4=l,i01,10 

CALL  EQUIV ( BUF ( L4 +  1 ),  1H+) 

160  CONTINUE 

180  GRIDP1=V1 (M) /FEXP 

JJ- (GRIDPl-GRIDLO) /DELTA+2 . 5 
IF(JJ-102)190,190,185 
185  JJ=102 

190  GRIDP2=V2 (M) /FEXP 

LL= (GRIDP2-GRIDLO) /DELTA+2 . 5 
IF(LL-102)195,195,192 
192  LL=102 

195  SAVF.1=BUF(JJ) 

SAVE2=BUF{LL) 

CALL  EQUIV(BUF(LL)  ,1H*) 

CALL  EQUIV { BUF ( JJ) ,1H=) 
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WRITE  (6,901)  GRIDP1,GRIDP2, (BUF(I),I=2,102) 

BUF(JJ)=SAVE1 
3UF(LL)=SAVE2 
200  CONTINUE 
C 

C  Scale  is  now  printed 
C 

TEND=DELTA*10.0 
BUF(l)=GRIDLO 
DO  220  1=2,11 
11=1-1 

BUF(I)=TEND+BUF(II) 

220  CONTINUE 

WRITE  (6,902)  (BUF(I) , 1=1, 11) 

C 

900  FORMAT (1H0,3X,1A6,5H  NO. , 13, 40X, 14HSCALE  FACTOR  =, 1PE9 . 1/lHO, 

+4X, 1H=, 6X, IH*) 

901  FORMATdH  ,  2F7 . 3,  IX,  102A1 ) 

902  FOPJ4AT(llX,  10  (3X,  F5 . 2, 2X)  ,  IX,  F5. 2) 

C 

RETURN 

END 

0======================================================================= 

C  FUNCTION  RDOT 

C - 

FUNCTION  RDOT  (A, MA, B, MB, N ) 

DIMENSION  A(MA, 1) ,B(MB, 1) 

RDOT  =0.0 
DO  10  J  =  1,N 

10  RDOT  =  RDOT  +  A(1,J)  *  B(1,J) 

RETURN 

END 

C  FUNCTION  CDOT 

C - 

FUNCTION  CDOT  (A, B,N) 

DIMENSION  A(l) ,B(1) 

CDOT  =0.0 
DO  10  I  =  1,N 

10  CDOT  =  CDOT  +  A(I)  *  B(I) 

RETURN 

END 

C  SUBROUTINE  CLIM 

C - 

SUBROUTINE  CLIM  (A, B, C, N, IRO) 

DIMENSION  A(l) ,B(1) 

DO  50  J  =  1,N 
IF  (J-IRO)  20,50,20 
20  A(J)  =  A(J)  -  B(J)  *  C 
50  CONTINUE 
RETURN 
END 

C  SUBROUTINE  SMOOTH 

C  For  each  iteration,  SMOOTH  calculates  new  window  widths  such  that 
C  the  final  statistical  error  (after  folding  with  the  window 
C  function)  should  converge  to  a  constant  specified  by  the  user. 
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C 


SUBROUTINE  SMOOTH ( ERROR, ERRF, WINW, NW, WIND, BLAB) 
DIMENSION  WINW(NW) ,ERRF(NW) , WIND (NW) , BLAB (NW) 
DIMENSION  WINDOW (130) 


C 

C(ENERGy)=40.*EXB(-.08*ENERGY)  +  50.*EXB(--  7*ENERGY) 

C 

C  For  each  iteration,  smooth  calculates  new  window  widths  such  that 
C  the  final  statistical  error  (after  folding  with  the  window 
C  function)  should  converge  to  a  constant  specified  by  the  user. 

C 

ERRM  =0.08 
DO  100  J=1,NW 

ERRF(J)  =ABS(ERRF(J) ) 

ENERGY=BLAB( J) *1.0E-6 
WMAX  =  1.5  *  WIND(J) 

IF ( ERRF ( J) .GT. ERRM)  GO  TO  10 
C 

WINDOW(J)=l/(  1/WINW(J)  + (1/C (ENERGY) ) *ALOG(ERROR/ERRF(J) ) ) 

GO  TO  20 
C 

10  WINDOW (J)  =  WMAX 

20  WINDOW(J)  -  ABS (WINDOW (J) ) 

100  I F (WINDOW ( J) .GT. WMAX)  WINDOW ( J) =WMAX 
C 

C  To  prevent  oscillations  of  the  window  widths,  the  results  are 
C  now  smoothed  with  a  5-point,  Ist-order  least  squares  fit  by  SE15. 
C 

CALL  SE15 (WINDOW, WINW,NW,IER) 

C 

RETURN 

END 

SUBROUTINE  SE15 


C 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 


PURPOSE 

TO  COMPUTE  A  VECTOR  OF  SMOOTHED  FUNCTION  VALUES  GIVEN  A 
VECTOR  OF  FUNCTION  VALUES  WHOSE  ENTRIES  CORRESPOND  TO 
EQUI DISTANTLY  SPACED  ARGUMENT  VALUES. 

USAGE 

CALL  SE15{Y,Z,NDIM,IER) 

DESCRIPTION  OF  PARAMETERS 

Y  -  GIVEN  VECTOR  OF  FUNCTION  VALUES  (DIMENSION  NDIM) 

2  -  RESULTING  VECTOR  OF  SMOOTHED  FUNCTION  VALUES 

(DIMENSION  NDIM) 

NDIM  -  DIMENSION  OF  VECTORS  Y  AND  Z 
lER  -  RESULTING  ERROR  PARAMETER 

lER  =  -1  -  NDIM  IS  LESS  THAN  5 

lER  =  0  -  NO  ERROR 

REMARKS 

(1)  IF  IER=-1  THERE  HAS  BEEN  NO  COMPUTATION. 

(2)  Z  CAN  HAVE  THE  SAME  STORAGE  ALLOCATION  AS  Y.  IF  Y  IS 
DISTINCT  FROM  Z,  THEN  IT  IS  NOT  DESTROYED. 

SUBROUTINE  AND  FUNCTION  SUBPROGRAMS  REQUIRED 
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NONE 


c 

C 

C  METHOD 

C  IF  X  IS  THE  (SUPPRESSED)  VECTOR  OF  ARGUMENT  VALUES,  THEN 

C  EXCEPT  AT  THE  POINTS  X (1 ) , X (2 ) , X (NDTM-1 )  AND  ” (NDIM) ,  EACH 

C  SMOOTHED  VALUE  Z(I)  IS  OBTAINED  BY  EVALUATING  AT  X(I)  THE 

C  LEAST- SQUARES  POLYNOMIAL  OF  DEGREE  1  RELEVANT  TO  THE  5 

C  SUCCESSIVE  POINTS  (X { I +K) , Y ( I+K) )  K  =  -2,-1,.,. ,2.  (SEE 

C  HILDEBRAND,  F.B.,  INTRODUCTION  TO  NUMERICAL  ANALYSIS, 

C  MC  GRAW-HILL,  NEW  YORK/TORONTO/LONDON,  1956,  PP.  295-302.) 

C - 

SUBROUTINE  SE15(Y, Z, NDIM, lER) 

DIMENSION  Y (NDIM) , Z (NDIM) 

Test  of  dimension 

IF(NDIM-5) 3,  1, 1 

Prepare  loop 

1  A=Y(1)+Y(1) 

C=Y(2) +Y(2) 

B=.2* (A+Y(l)+C+Y(3)-Y(5) ) 

C=. 1* (A+A+C4Y(2) +Y(3) +Y(3) +Y(4) ) 

Start  loop 

DO  2  1=5, NDIM 
A=B 
B=C 

C=.2* (Y(I-4)+Y(I-3)+y(I-2)+Y(I-l)+y(I) ) 

2  Z(I-4)=A 
C 

C  End  of  loop 

C 

C  Update  last  four  components 

C 

A=Y(NDIM) +Y(NDIM) 

A=.l* (A+A+Y(NDIM-1) fY (NDIM-1 ) +Y (NDIM-1 ) +y (NDIM-2 ) +Y (NDIM-2 ) 

+  +Y(NDIM-3)) 

Z(NDIM-3)=B 
Z (NDIM-2 )=C 
Z (NDIM-1) =A 
Z (NDIM)=A+A-C 
IER=0 
RETURN 
C 

C  Error  exit  in  case  NDIM  is  less  than  5 
C 

3  IER=-1 
RETURN 
END 

C  end  of  FORIST.FOR 


Appendix  F:  GRAPHICS. TK  Computer  Model 


.  '  VAHiXBLe  SHEET  . . .  ■  for  Acadeaic  Us*  Only 

St  Input -  Noae - -  Output -  Unit -  Couasnt - 

Graphics. TX  23  Oct  93 

by  Robert  S.  Pop* 

A  PROGRAM  TO  GRAPH  THE  RESULTS  FROM 
THE  rORIST  UHFOLOING  CODE 


energy  -  input  fro*  TOGRAPH.OUT  (*7) 
•nargy  -  converted  for  x~axls  (M*V) 
spectrua,  unsaoothed 
error,  unsaoothed 

high  bound  on  spectrun  error,  unssooth 
low  bound  on  spectrua  error,  unsBooth 

energy  -  Input  fros  TOGRAPH.OUT  (*V) 
energy  -  converted  for  x-axis  (MeV) 
spectrua,  saoothed 
error,  smoothed 

high  bound  on  spectrua  error,  saooth 
low  bound  on  spectrua  error,  saooth 


. .  I  RULE  SHEET  n  ■  ..i— .n.  for  Acadeaic  Us*  Only 

S  Rule - — - - - 

"Convert  energy  froe  eV  to  MeV 

•  *1  -  El  /  le6 

•  *0  -  EO  /  1*< 

"Coaput*  error  bounds  on  spectrua 

•  SpecUpl  >  UTBl  *  EPSLOHl 

•  SpecOnl  •  UTBl  •  EPSLOHl 

•  SpecUpO  -  UTBO  *  EPSLOHO 

•  SpecDnO  •  UTBO  -  EPSLOHO 


II  LIST  SHEET  '  I  I  — .  for  Acadeaic  Use  Only 

Naae -  Eleaents —  Unit -  Coaaent - 

I.  Us*  "/SAL*  (ASCII  load)  coaaand  to  Input 
lists  froa  TOGRAPH.OUT 
3.  Hit  PIO  to  create  reaaining  lists 
3.  Plot  result* 


L 

L 

L 

L 

L 

L 

L 

L 

L 

L 

L 

L 


El 

*1  0 

UTBl 
EPSLOHl 
SpecUpl  0 
SpecDnl  0 

BO 

eO  0 

UTBO 
EPSLOHO 
SpecUpO  0 
SpecOnO  0 


Naae - 

unsaoothl 
saooth 0 


. .  ■—  PLOT  SHEET  ' 

Plot  Type —  Display  Option-  Output  Device- 
Line  chart  l.EGA 
Line  chart  I . EGA 


>  for  Acadeaic  Us*  Only 

Title - 

Add  title 
Add  title 
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Appendix  G:  ZSHIFT.OIJT  Data  File 


This  is  the  file  ZSHIFT.OUT  created  by  ZSHIFT.  This  file  contains  each  of  the 
four  original  spectra,  the  four  corresponding  zero-shifted  spectra,  and  the  four  final  200- 
channel  zero-shifted  spectra  for  input  to  FORIST.  The  purpose  of  this  file  is  to  serve  as  a 
debugging  tool  for  ZSHIFT.  The  data  in  this  file  corresponds  to  the  other  computer 
examples  in  the  appendices.  This  file  is  presented  in  lanscape  format  in  order  to  preserve 
the  form  of  the  actual  output  from  ZSHIFT,  which  is  formatted  for  output  on  a  wide- 
carriage  line  printer. 
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NAME  OF  RESPONSE  MATRIX  USED  IN  PROBLEM  =UINEUT.MAT 
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Appendix  H:  SPECTRA  IN  Data  File 


The  following  is  the  data  file  called  SPECTRA.IN,  which  is  created  by  ZSHIFT  to 
be  the  input  file  for  FORIST.  This  version  of  the  file  corresponds  to  the  other  computer 
files  contained  in  the  appendices.  The  format  of  this  file  is  important  and  must  be 
maintained  if  FORIST  is  to  be  able  to  read  the  file  properly.  The  page  breaks  of  this  file 
have  been  included  so  that  the  output  would  fit  into  the  format  of  the  thesis.  In  the  acutal 
file,  there  are  no  blank  lines  between  the  pages. 
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Appendix  I:  FQRIST.OUT  Data  File 


This  appendix  contains  the  output  file  from  FORIST,  called  FORIST.OUT,  as  well 
as  a  key  to  the  information  that  is  contained  on  each  page  of  this  output  file.  Because  the 
output  file  is  meant  to  be  printed  on  a  wide-carriage  printer,  FORIST.OUT  is  presented  in 
ianscape  format.  To  correspond  with  the  key,  page  number  headings  have  been  added  to 
indicate  the  beginning  of  each  page.  The  page-by-page  key  to  this  output  file,  which  is 
presented  in  Table  14,  has  been  adapted  from  Johnson. 

[15:16-18] 


Table  14.  Description  of  contents  of  FORIST.OUT  file 


Pace  Number 

Description 

1 

Title  or  response  matrix,  FORIST  input  parameters 

2 

Binning  table 

3 

Normalizations,  gains,  and  nominal  gains 

4 

Foreground  data,  high  and  low  gain 

5 

Rearranged  high-gain  foreground  data 

6 

Rearranged  low-gain  foreground  data 

7 

Background  data,  high  and  low  gain 

8 

Rearranged  high-gain  background  data 

9 

Rearranged  low-gain  background  data 

10 

Rebinned  high-gain  foreground  data 

11 

Rebinned  low-gain  foreground  data 

12 

Rebinned  high-gain  background  data 

13 

Rebinned  low-gain  background  data 

14-15 

Summary  of  binned  input: 

COUNTS:  grouped  counts/light  unit/second 

ERROR:  statistical  error  in  COUNTS 

LOCH  AN:  lower  channel  of  bin 

UPCHAN:  upper  channel  of  bin 

LO  P.  H.:  lower  pulse  height  boundary  of  bin 

UP  P.  H.:  upper  pulse  height  boundary  of  bin 

F.  A.F. :  Forte  Acceptance  Factor 

Table  14,  (continued) 


Paste  Number 

Description 

16-17 

Overlap  parameters  SCOFF,  FUDGE,  and  REFLU 

Plot  of  binned  counts  and  errors 

18-19 

Binned  input  adjusted  for  smooth  overlap;  chi-square  test 

I:  bin  number 

VLO:  lower  pulse  height  boundary  of  bin 

B;  COUNTS  adjusted  for  smooth  overlap 

S;  ERROR  adjusted  for  smooth  overlap 

BADJ:  refolded  vector 

PCTDEV:  relative  error,  in  percent 

ESD:  estimated  standard  deviation 

20 

Indicators  foi  va!id  problem 

1 :  matrix  i:;v’''sion  accuracy 

2:  consistency  with  non-negativity 

3:  chi-square  consistency 

21-22 

Unfolded  solution  summary  ("zeroth"  smoothing  iteration) 

ENERGY:  energies,  in  eV 

PLO:  lower  bound  (-1  std.  dev.)  of  unfolded  spectrum 

PUP:  upper  bound  (+1  std  dev.)  of  unfolded  spectrum 

PAVE:  unfolded  spectrum,  in  neutrons/cm^2-sec-MeV 

PCT  W:  window  width,  in  percent 

ERRl :  statistical  error  in  PAVE 

ERR2:  numerical  error  in  PAVE 

23-24 

Pseudo-semilog  plot  of  PLO  and  PUP  for  "0th"  smoothing  iteration 

25-26 

Unfolded  solution  summary  (first  smoothing  iteration) 

27-28 

Pseudo-semilog  plot  of  PLO  and  PUP  for  first  smoothing  iteration 

29-30 

Linear  plot  of  PLO  and  PUP  for  final  smoothing  iteration 

31-32 

Summary  of  unfolded  solution  before  smoothing 

ELAB:  energies,  in  eV 

X:  solution  vector 

Q:  weak  upper  bound  for  X 

ERRJ 1 :  statistical  error  in  X 

UHAQJ:  numerical  error  in  X 

The  file  FORI  ST.  OUT  corresponding  to  the  rest  of  the  code  examples  in  the 
apendices  will  now  be  presented  in  its  entirety. 
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J:  TOGRAPH  OUT  Data  File 


The  file  TOGRAPH. OUT  is  generated  by  FORIST  and  is  used  as  the  data  input  to 
the  graphing  routine  GRAPHICS  . TK.  This  fde  contains,  for  each  i*eration  of  FORIST, 
the  unfolded  spectrum  (UTBx),  the  energy  scale  of  the  unf  jlded  spectrum  (Ex),  and  the 
uncertainty  in  the  unfolded  spectrum  (EPSLONx).  The  "x"  in  each  variable  name  is  an 
integer  which  begins  with  its  largest  value  corresponding  to  the  unsmoothed  spectrum  arid 
decreases  by  one  with  each  smoothing  iteration  until  zero  is  reached.  If  FORIST  is  run  in 
its  ususal  mode  with  one  smoothing  iteration  (ITER=1),  there  will  be  a  total  of  six  data 
lists  in  this  file.  The  lists  in  this  file  correspond  to  the  computer  examples  in  the  other 
appendices. 

This  data  file  is  created  by  FORIST  to  oe  in  the  format  required  by  TK  Solver 
Plus.  Each  data  set  begins  with  the  name  of  the  variable,  which  is  followed  by  a  colon. 

The  data  is  then  listed  one  point  at  a  time  with  a  comma  following  each  data  point  but  the 
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procedure  are  on  the  order  of  5-10%.This  thesis  is  intended  to  be  a  summary  of  the  relevant  literature 
and  a  user’s  guide  to  the  PSD  spectrometer. 
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